A probabilistic assessment of calcium carbonate export and dissolution in the modern ocean

Abstract. The marine cycle of calcium carbonate (CaCO3) is an important element of the carbon cycle and co-governs the distribution of carbon and alkalinity within the ocean. However, CaCO3 export fluxes and mechanisms governing CaCO3 dissolution are highly uncertain. We present an observationally-constrained, probabilistic assessment of the global and regional CaCO3 budgets. Parameters governing pelagic CaCO3 export fluxes and dissolution rates are sampled using a Monte 5 Carlo scheme to construct a 1000 member ensemble with the Bern3D ocean model. Ensemble results are constrained by comparing simulated and observation-based fields of excess dissolved calcium carbonate (TA∗). The minerals calcite and aragonite are modelled explicitly and ocean–sediment fluxes are considered. For local dissolution rates either a strong or a weak dependency on CaCO3 saturation is assumed. In addition, there is the option to have saturation-independent dissolution 10 above the saturation horizon. The median (and 68 % confidence interval) of the constrained model ensemble for global biogenic CaCO3 export is 0.90 (0.72–1.05) Gt-C yr−1, that is within the lower half of previously published estimates (0.4–1.8 Gt-C yr−1). The spatial pattern of CaCO3 export is broadly consistent with earlier assessments. Export is large in the Southern Ocean, the tropical Indo–Pacific, the northern Pacific and relatively small in the Atlantic. The constrained results are 15 robust across a range of diapycnal mixing coefficients and, thus, ocean circulation strengths. Modelled ocean circulation and transport time scales for the different setups were further evaluated with CFC11 and radiocarbon observations. Parameters and mechanisms governing dissolution are hardly constrained by either the TA∗ data or the current compilation of CaCO3 flux measurements such that model realisations with and without saturation-dependent dissolution achieve skill. We suggest 20 to apply saturation-independent dissolution rates in Earth System Models to minimise computational costs.

typically neglected in previous studies on open water CaCO 3 dissolution, is arising from oceansediment interactions and the influence of associated burial and redissolution fluxes on alkalinity and carbon concentrations in the open ocean (Archer, 1996).
Here we propose an alternative, probabilistic assessment of the global CaCO 3 budget. We consider explicitly the transport and mixing of alkalinity and account for ocean-sediment interactions to probabilistically constrain the CaCO 3 cycle with the observation-based distribution of the TA * tracer. TA * reflects the imprint of the CaCO 3 cycle on alkalinity (see Sect. 2) and is therefore im- concentrations of TA * . The most probable export and dissolution fluxes are then the ones resulting in modelled TA * fields that match closely the observation-derived TA * distribution. The aim is to assign uncertainty estimates that will be consistent with both the observations and model equations -a classical data assimilation problem (see Sect. 3.3). It is expensive to perform simulations with interactive sediments (Heinze et al., 1999;Gehlen et al., 2006;Tschumi et al., 2011) due to long spin-80 up times required to bring ocean and sediments in equilibrium and we account a posteriori for the influence of sediment burial and dissolution fluxes on TA * (see Sect. 3.3). This approach represents an alternative to the interpretation of concentrations without a physical model or the interpretation of flux/rate measurements within the water column, which are generally much sparser (∼ 156 flux measurements in the water column (Wilson et al., 2012) and ∼ 56 benthic dissolution flux measurements 85 (Berelson et al., 2007), globally) and more difficult to obtain. We apply the Wilson et al. (2012) data compilation of water column fluxes as an additional constraint for comparison. Additional sensitivity analyses with respect to vertical diffusion (k dia , see Sect. 4) are illustrated. The results are compared and contrasted to data-based estimates of export and dissolution (as summarised in Berelson et al., 2007). Finally, implications for the parameterisation of CaCO 3 dissolution in Earth System Models 90 are discussed.
2 Observation-derived TA * Total alkalinity data are from the GLODAP carbon climatology (Key et al., 2004) and salinity (Antonov et al., 2010), temperature (Locarnini et al., 2010), oxygen and phosphate (Garcia et al., 2010a, b) data are from the World Ocean Atlas. These gridded data-products represent objectively 95 analyzed climatological fields of the respective oceanic variables and are based on samples taken during the previous decades. They serve to split the alkalinity signal into its different physical and biogeochemical components such as TA * , our target variable in the data assimilation. We first re-gridded all required, gridded datasets to the Bern3D model grid (40x41x32 grid boxes) using the area-weighted regridding method of Ferret before deriving the other properties. 100 TA * , our target variable in the data assimilation, is a constructed tracer (Feely et al., 2002;Sabine et al., 2002a;Chung et al., 2003;Koeve et al., 2014) to exclusively capture the imprint of CaCO 3 dissolution on alkalinity. TA * is, in this sense, one of three components of measured total alkalinity (TA, Eq. 1, mean concentration 2.427 mol m −3 based on the regridded GLODAP dataset). TA * can be extracted from TA by accounting for preformed (TA 0 , Eq. 2) and remineralised alkalinity (TA r ,105 Eq. 3): TA = TA 0 + TA r + TA * . (1) TA 0 is the background or preformed concentration, set at the ocean surface, and mixed conservatively in the ocean interior. Accordingly, TA r and TA * are by definition zero in the surface ocean.
TA 0 = a 0 + a 1 · S + a 2 · T + a 3 · PO (2) 115 The coefficients, a i , are estimated from observations of TA (Key et al., 2004), S (Antonov et al., 2010), T (Locarnini et al., 2010), and PO (Garcia et al., 2010a, b) in the surface ocean. The linear regression fit is sometimes further subdivided to include only specific basins (Feely et al., 2002;Koeve et al., 2014). Table 1 summarises different regressions, including previously published regression estimates as well as new estimates by us calculated on the Bern3D model grid (B3D). The 120 root mean squared deviations (RMSE) between these fits are smaller than 0.00465 mol m −3 (i.e. smaller than 0.2 %, excluding the Gruber et al. (1996) equation, which relies on an older database from 1996). From these different options, we accordingly chose the Global B3D linear regression including 3 explanatory variables and global datasets as robust regression for TA 0 which yields a mean concentration of ∼ 2.389 mol m −3 .

125
TA r is linked stoichiometrically to the apparent oxygen utilisation (AOU, Garcia et al., 2010b, Eq. 3) and accounts for decreases in TA due to the oxidation of organic nitrogen, phosphorous, and sulfur (OM).
In addition, we note that AOU has been suggested to overestimate true oxygen utilisation by 20-25 % (Ito et al., 2004;Duteil et al., 2013). Accordingly, TA r might be associated with an uncertainty 135 of ∼ 20 %.
The remaining signal, then, is TA * (mean concentration ∼ 57 mmol m −3 based on GLODAP and our reference choices to derive TA 0 and TA r ), the changes in alkalinity due only to the CaCO 3 cycle ( Fig. 1). The global average RMSE of any of the described ways (16 in total) of accounting for TA r (25 % lower AOU or different r Alk:OM ) and TA 0 (either 2 or 3 explanatory variables, and either global 140 or regional) from our reference choices is ∼ 4 mmol m −3 (average RMSE 3.9, 3.8, 3.1 mmol m −3 in the Atlantic, Pacific and Indian Ocean, respectively, i.e. ∼ 7 %). It is to note that this approach, inherent to its empirical nature, yields slightly negative TA * values in some places. TA * integrates to about 37.5 Pmol C or 75 Pmol Alk-equivalents of which ∼ 41 % come to lie above the calcite saturation horizon (similar to Koeve et al. (2014) who find 44.7 % of the TA * inventory above the 145 calcite saturation horizon). These estimates are robust across the different sources of uncertainty.
TA * concentrations are expressed in alkalinity equivalents throughout this paper and we do not divide TA * concentrations by a factor of 2 as done in many observational studies which express TA * in terms of carbon changes. TA * inventories are in Pmol C and CaCO 3 fluxes are given in carbon units (Gt-C yr −1 , mmol-C m −3 yr −1 , mmol-C m −2 yr −1 ).

The modelling framework
To constrain the alkalinity fluxes associated with CaCO 3 cycling, we introduce different export and

The Bern3D Model
The Bern3D model couples a dynamic ocean, sea ice, an energy-moisture balance atmosphere, a marine biogeochemical cycle, a dynamic global vegetation model, and an ocean sediment module (Müller et al., 2006;Tschumi et al., 2011;Ritz et al., 2011;Roth et al., 2014). Here an ocean version with a horizontal resolution of 41 by 40 grid cells and 32 logarithmically-scaled vertical 170 layers is used (see also Roth et al., 2014). The horizontal resolution is the same for the components atmosphere, ocean, sea ice, and sediments of the Bern3D model. Transport and mixing of tracers in the ocean is based on Edwards et al. (1998) and Müller et al. (2006) as a three-dimensional frictional geostrophic model. The model has an isopycnal diffusion scheme and Gent-McWilliams parameterisation for eddy-induced transport (Griffies, 1998). The NCEP/NCAR monthly wind-stress 175 climatology (Kalnay et al., 1996) is prescribed at the surface. Air-sea gas exchange for CO 2 is implemented according to OCMIP-2 protocols (Najjar et al., 1999;Orr and Najjar, 1999). The global mean air-sea transfer rate is reduced by 19% compared to OCMIP-2 to match observation-based estimates of natural and bomb-produced radiocarbon . A two-dimensional energy moisture balance model represents the atmosphere (Ritz et al., 2011 Parekh et al. 190 (2008) and Tschumi et al. (2011). One third (33 %) of the new production is exported out of the euphotic zone (defined at 75 m) as particulate organic matter (POM), with the reminder contributing to the dissolved organic matter pool. Within biogeochemically similar regions/provinces a parameter termed rain ratio linearly scales the pattern of POM export to CaCO 3 export (silica-limitation is not considered here). We define eight such regions, each assigned an independent value for the In the model, we implemented TA * as an explicit, idealised tracer (Koeve et al., 2014). It captures the alkalinity equivalents of CaCO 3 dissolution whenever it occurs and mixes this signal, accordingly. TA * values are set to zero throughout the surface ocean.

205
In sensitivity simulations with the sediment module enabled, the flux of CaCO 3 , and of other particles, reaching the seafloor is passed to the sediment module from where a fraction potentially re-dissolves back into the water column. In simulations without the sediment module, the entire flux reaching the ocean floor re-dissolves back into the water column. Simulated TA * concentrations tend to be lower with the sediment module enabled than without sediment module, because a fraction of 210 the CaCO 3 export flux is removed from the ocean and buried in the geosphere. The sediment diagenesis model (Heinze et al., 1999;Gehlen et al., 2006;Tschumi et al., 2011) features the same horizontal resolution as the ocean model and 10 layers resolving the top 10 cm of the seafloor. It dynamically calculates the transport, remineralisation/redissolution and bioturbation of solid material within the top 10 cm of the seafloor as well as pore-water chemistry and diffusion as described in 215 detail in Tschumi et al. (2011). Solutes diffuse over a boundary layer of 1 cm between the sediment column and the lowermost ocean grid cell.

CaCO 3 dissolution within the water column
The dissolution equations for calcite and aragonite are implemented using equations with identical 235 functional forms, but with different parameters for each mineral. In the following, we do not explicitly distinguish the two minerals to ease notation. The dissolution of calcite and aragonite below the euphotic zone is assumed to be a function of the saturation state of the bulk seawater, Ω, and the par-ticle concentration per unit of water volume, [CaCO 3 ] (see Gangstø et al. (2011) for a discussion): 240 k eff denotes a first order rate constant for either calcite or aragonite. It is defined as: k 0 and k bg are rates in units of 1/time. H is the Heaviside function, which is zero for supersaturated and 1 for undersaturated water. The saturation state is defined by the ratio of the product of calcium ion concentration times carbonate ion concentrations to the saturation product, K sp (Mucci, 245 1983): Ω values larger than one correspond to oversaturated and Ω values smaller than one to undersaturated conditions. Ω and thus k eff are grid cell specific. At super-saturation, the dissolution rate k eff equals the constant background rate k bg (which can be zero). With increasing under-saturation, the dissolu-250 tion rate increases towards its maximum value (k bg + k 0 ). n is a unitless parameter and determines the deviation from linearity of this increase. For simplicity and to avoid the addition of further free parameters, a constant sinking velocity, v, is assumed and identical for both calcite and aragonite particles. The flux profile of CaCO 3 then takes the form where F is the downward flux of either calcite or aragonite particles per unit area evaluated at the bottom of each tracer grid cell at depth z k . i, j, and k are grid cell indices indicating longitude, latitude, and depth. ∆z k denotes grid cell height. The export flux is set equal to F (z = 75 m) at the depth of the euphotic zone. Particles are dissolved instantaneously and sinking is not explicitly resolved in this formulation, reducing computational costs. v/k eff is in units of m, and, if assumed constant, can 260 be interpreted as the dissolution length scale, i.e. the depth at which the flux has decreased to 1/e (to ∼ 37 %) of the export flux at 75 m. Any particles reaching the sea floor are dissolved completely in the appropriate lowermost box, except when the sediment module is included.
In Fig. 2 we illustrate three dissolution cases to explore the sensitivity of TA * to dissolution profiles which cover the sampled uncertainty in dissolution rate parameters. only below the saturation horizon. In the fast (slow) case, aragonite (chosen to represent 10 % of the total) and calcite (90 % of total export) dissolves quickly (slowly) below the saturation horizon.

270
Aragonite and calcite dissolve within a few hundred meters below the saturation horizon in the fast case, while most CaCO 3 dissolves on the ocean floor in the slow case.
The choice of the dissolution rate profile has a substantial influence on the simulated TA * inventory ( Fig. 2, Table 2). The global TA * inventory is 38 Pmol C for the case with a constant, saturationindependent dissolution and a rain ratio of ∼ 7 %, which is close to the observation-derived inventory 275 of 37 Pmol C. The simulated inventory is 48 and 63 Pmol C for the fast and slow cases (where no dissolution occurs above the saturation horizon), respectively, and by that substantially higher. TA * accumulates in the deep ocean when all CaCO 3 is dissolved below the saturation horizon of aragonite and calcite and no dissolution is permitted above. The TA * inventory and concentrations are sensitive to the choice of the dissolution rate profile, supporting our choice of TA * as a target variable 280 to constrain dissolution rates.

Skill scores
Global skill scores, S m , are assigned to each member, m, of the Latin-Hypercube ensemble:

320
MSE rel is the relative mean squared error of the simulated TA * concentrations from member m with respect to observation-derived TA * : The sum includes all grid cells (indexed j). TA * model The above explanations apply for any simulated quantity of interest. In the following we will present PDFs, median values, and 1σ confidence ranges for aragonite and calcite export and dissolution as well as for tracer concentrations at individual grid cells or integrated over regions or 355 the whole ocean. Spatial integrations are done for each ensemble member individually and before computing the PDFs and associated measures from the full ensemble.

A first order correction for ocean-sediment TA * fluxes
CaCO 3 burial removes alkalinity from the ocean water column and lowers concentrations and the overall TA * inventory relative to a run without the sediment module. Riverine input compensates this 360 loss. This input of alkalinity is added to the surface ocean and by that to the part of the preformed alkalinity component (TA 0 ), leaving TA * unchanged. CaCO 3 export and dissolution within the water column and the corresponding fluxes of alkalinity and TA * remain largely unchanged between runs with and without sediments.
Ideally, the ensemble would be run fully interactively with the sediment module enabled to ac-365 count for all important processes within the CaCO 3 cycle. However, this is computationally too expensive as the sediment module requires a long spin up to achieve equilibrium. A first order correction term TA * sedcorr that accounts for the influence of CaCO 3 burial and dissolution fluxes on TA * is estimated as follows. First, the skill scores are computed as described above with TA * sedcorr set to zero. Then, the 14 best ensemble members are selected and re-run with the sediment module runs without the sediment (see Table 3: 0.9 Pmol C in the Atlantic, 4.6 Pmol C in the Pacific, and 1.4 Pmol C in the Indian Ocean). This reduction is equal to the inventory of the TA * sedcorr correction and ∼ 16 and ∼ 20 % of the observation-derived TA * inventory of the Pacific and Indian Ocean, 380 respectively. TA * concentrations are affected relatively uniformly below 1 km and differences in TA * between simulations with and without sediment module tend to vanish toward the surface ocean ( Figs. 3 and 4). Correspondingly, the spatial patterns of TA * are very similar for simulations with and without sediments; correlation coefficients are > 0.99. We expect that the correction will tend to increase the export flux of CaCO 3 in the optimisation to compensate for the loss by burial, but 385 will not strongly affect dissolution parameters as the spatial patterns of TA * remain similar with or without correction.
A few caveats apply to this first order estimate of TA * sedcorr . First, our approach involves a few, likely minor, technical inconsistencies. Aragonite is not treated explicitly in the sediment module and all CaCO 3 is assumed to be in the form of calcite; this may tend to bias sediment burial high 390 as calcite is less soluble. The saturation state within the sediments is computed interactively with modelled ocean boundary conditions; this may locally lead to inconsistencies as CaCO 3 dissolution within the ocean water column is computed using the prescribed, observation-based saturation state.
The alkalinity flux associated with organic matter remineralisation within the sediment is not explicitly distinguished and included in the flux of TA * from sediments to the ocean; this results in a bias 395 on order of 5 % in the redissolution flux and a negligible influence on the sediment correction. could potentially be reduced by introducing more than one background dissolution rate constant (k bg ) or a depth-dependent particle sinking velocity. However, this may simply mask deficiencies in the circulation and we do not attempt such a solution. Generally, the correlation between observed and modelled TA * is remarkably high.

Probabilistic estimates of CaCO 3 and alkalinity fluxes 435
The estimated global median export flux of CaCO 3 at 75 m -with its 68 % confidence interval -is 0.90 (0.72-1.05) Gt-C yr −1 (Table 4 and Fig. 6). Basin-wide, we find CaCO 3 median export fluxes (and 68 % confidence intervals) of 0.12 (0.078-0.17) Gt-C yr −1 from the Atlantic, 0.55 (0.39-0.68) Gt-C yr −1 from the Pacific, and 0.23 (0.14-0.32) Gt-C yr −1 from the Indian Ocean (Table 4 and Fig. 6). Regionally, largest export fluxes are simulated in the Southern Ocean sector of the 440 Pacific, in the Pacific equatorial upwelling regions, and in the northwestern Pacific (Fig. 7). In the Indian Ocean export fluxes are highest in its eastern tropical regions and in its section of the Southern Ocean. In the tropical and northern Atlantic export fluxes are generally low, consistent with the low TA * values in the bulk of NADW and AABW.
Median CaCO 3 export production per unit area in the constrained ensemble is considerably 445 lower in the Atlantic sector of the Southern Ocean (135 mmol-C m −2 yr −1 ) as compared to the Pacific (301 mmol-C m −2 yr −1 ) and Indian (326 mmol-C m −2 yr −1 ) sectors (Fig. 7 top). This is attributable to the choice of the regional boundaries for the rain ratio regions and the assumption that the spatial pattern of export within a region is identical to the pattern simulated by the standard version of the model. The standard model yields relatively little zonal variation in the CaCO 3 export 450 fluxes in the Southern Ocean in contrast to the data assimilation with lower than zonally-averaged export in the Atlantic. This reflects the much lower TA * reconstructed in the deep Atlantic as compared to the deep Pacific and Indian (Fig. 1). A large export in the Atlantic sector of the Southern Ocean tends to yield high simulated TA * concentrations in the Antarctic Bottom Water that fills the deep Atlantic. The Monte Carlo data assimilation therefore requires low CaCO 3 export in the Atlantic 455 sector to minimize model-data mismatches in the deep Atlantic. It is difficult to correctly represent water mass formation and circulation in the Southern ocean and our model may be biased. A known bias is that the Atlantic Bottom Water circulation is too sluggish, also evidenced by simulated low radiocarbon signatures (Figs A2. A4). The influence of a potential bias in South Atlantic export on global CaCO 3 export is estimated to be relatively small; assuming the same (median) CaCO 3 ex-460 port per unit area in the Atlantic sector as estimated for the Indian sector would yield 0.06 Gt-C yr −1 higher export than suggested by the ensemble median. While our Monte Carlo approach is suitable to estimate export fluxes over larger regions, the detailed spatial patterns in CaCO 3 export remain unconstrained.
Average dissolution profiles for aragonite (red) and calcite (olive) in different ocean regions are displayed in Fig. 8. A peak in aragonite and calcite dissolution is located at or below the depth of the 475 aragonite and calcite saturation horizon, respectively. These dissolution peaks are associated with the saturation-dependent dissolution rate coefficients (Eq. 5). As is the saturation horizon, these peaks are located deep down in the water column of the north and tropical Atlantic region (at 3-4 km depth for aragonite), at intermediate depth of the South Atlantic, Indian, and Pacific (1.5 km for aragonite) and at relatively shallow depth of the tropical and north Pacific region (∼ 700 m for aragonite).

480
As the calcite saturation horizon is found deeper in the water column than the aragonite saturation horizon, so are these dissolution peaks located deeper down in the water column for calcite than for aragonite.
Our constrained ensemble includes non-zero values for the background dissolution rate. Consequently, calcite and aragonite are dissolving throughout the water column, irrespective of the sat-485 uration state in the Atlantic, Pacific and Indian Ocean. The percentage of the export flux which dissolves in waters supersaturated with respect to calcite are 72 % (60-80 %), 43 % (30-53 %), 68 % (55-78 %) in the Atlantic, Pacific and Indian Ocean, respectively. We will further investigate in Sect. 4.4 to which extent the finding that a fraction of the CaCO 3 export dissolves above the calcite saturation horizon and our export and dissolution flux estimates are robust. we set the sediment correction (TA * sedcorr ) of the TA * field to zero to calculate potential skill scores (Eq. 9). This alternative case also illustrates the potential error due to the neglect of ocean-sediment fluxes. Export fluxes of CaCO 3 are 2, 12, 8, and 8 % smaller in the Atlantic, Pacific, Indian and global ocean, respectively while the variance is about the same. This is not surprising as we have already noted (Sect. 3.3) that a certain burial flux tends to decrease the TA * pattern uniformly, i.e. 500 the simulations with and without the sediment correlate highly in terms of TA * . The PDF of fluxes is therefore shifted to higher export fluxes, while the preference for dissolution parameters remains the same. The CaCO 3 that dissolves above 1500 m is 4 % higher, 5 % lower, the same, and 4 % lower in the Atlantic, Pacific, Indian and global ocean, respectively, in the setup that neglects sediment fluxes compared to the standard setup. The reference simulation achieves a mean (natural) ∆ 14 C of DIC of −160 ‰. Correspondingly, ocean radiocarbon signatures become too low (−176 ‰) with the low diapycnal mixing rate and too high with the high (−126 ‰) diapycnal mixing rate (see Fig. A4). Simulated surface-to-deep ∆ 14 C gradients are too low (high) relative to the observed gradients for the high (low) diapycnal diffusivity parameter, thereby indicating too fast (slow) surface-to-deep water exchange (Fig. A4).
The global observation-based CFC11 inventory is estimated to 575 Mmol. The reference simulation yields 513 Mmol CFC11 as the mean inventory over the modelled period between 1990 and 2000. In the low mixing simulation, this inventory is even lower (479 Mmol) and in the high mixing 525 simulation it is too high (631 Mmol). On a global scale, our reference choice is therefore in better agreement with these physical tracers (see Fig. A4).
We vary k dia for three illustrative dissolution rate profiles introduced in Sect. 3.2 (see Fig. 2 and Table 2). In the low diapycnal diffusivity case, CaCO 3 export is 11 % lower and in the high diapycnal diffusivity case, CaCO 3 export is 30 % higher compared to the standard case. Larger Perhaps somewhat surprisingly, the simulated TA * inventory is relatively weakly affected by the choice of k dia ; the global TA * inventory varies by less than 7 % across the range of k dia (Table 2).
In other words, variations in the magnitude of ocean ventilation do hardly affect the TA * inventory for a given dissolution rate profile. Higher (lower) export under higher (lower) mixing compensate 540 each other. These changes in TA * inventory are smaller than the influence of the sediment correction (see Table 3) or the choice of the dissolution profile. In conclusion, simulated TA * is only weakly affected by uncertainties in the diapycnal mixing coefficient.
The choice of the dissolution rate profile has a substantial influence on the simulated TA * inventory ( Table 2). As mentioned previously, TA * accumulates in the deep ocean when all CaCO 3 545 is dissolved below the saturation horizon of aragonite and calcite and no dissolution is permitted above. This raises the question whether surface-to-deep transport is too slow in our model. As mentioned above, the radiocarbon signatures as well as CFC11 concentrations are on average close to observations (Figs. A3 and A2). Increasing ocean ventilation by increasing k dia results in too young radiocarbon signatures and higher than observed CFC11 concentrations. However, it does not sub-550 stantially reduce the overestimation of the TA * inventory in those cases.

How to parameterise CaCO 3 dissolution in an Earth System Model?
An important question is how to formulate the dissolution rate of calcite and aragonite particles in Earth System Models. Should the dissolution rate be a function of the simulated aragonite and calcite saturation state of the surrounding water? Should dissolution above the saturation horizon be 555 permitted? We analyse the relationship between model skill, dissolution parameterisation, CaCO 3 export, and TA * inventories in the model ensemble (see Fig. 9) to address these questions. To achieve a high skill (green-to-red coloured dots in the upper panels of Fig. 9), an individual ensemble member needs to reproduce the observation-based TA * inventory (dashed line) within a limited range.
We identify an export range within which TA * can be reproduced skillfully (vertical green range 560 in the upper panels of Fig. 9). Surprisingly, a high skill is achieved across the range of different dissolution schemes applied and for a broad range of parameter values. In other words, neither the dissolution scheme nor its parameters are well constrained by the observation-based TA * field. This is illustrated by plotting the value of the background dissolution rate, k bg (lower panels of Fig. 9) as a function of the TA * inventory and export. Generally, the higher the global CaCO 3 export flux, 565 the higher the background dissolution rate required to achieve a high skill. Likewise, lower export can be distributed skillfully without dissolution in supersaturated waters (k bg = 0). Apparently, there are trade-offs between the magnitude of export and the applied dissolution parameterisation in terms of TA * , suggesting that export and dissolution parameters can only be constrained simultaneously within limits when using TA * as the only constraint.

570
These findings are in line with the results from our sensitivity simulations (Fig. 2, Table 2). The parameterisation with dissolution above the saturation horizon (constant) yields the lowest TA * inventory, followed by the scheme with fast dissolution below the saturation horizon, and the scheme with dissolution near the ocean floor (slow). A high (low) export is thus required for parameterisations with high (low) dissolution above the saturation horizon to simulate the observation-based 575 TA * . Further, the different dissolution schemes yield highly correlated TA * fields; the correlation coefficient between the global fields from the fast and slow dissolution scheme is 0.84, and between the fields from the fast and constant scheme is 0.88 in simulations with the same CaCO 3 export (Fig. 2).
The average RMSE between the fast and slow dissolution scheme is 43.2 mmol m −3 and between the fast and constant is 24.5 mmol m −3 . The high spatial correlation in simulated TA * and uncer-580 tainties in CaCO 3 export make it difficult to distinguish different dissolution parameterisations. The magnitude of CaCO 3 exports modulates absolute TA * concentrations and thus model-data bias and root mean square errors. Given these uncertainties, we cannot objectively determine the preferred dissolution scheme. The TA * -CFC age method is heavily criticised by Friis et al. (2006) and tends to bias estimates systematically towards high values (Friis et al., 2006). This method and its shortcomings are further We find peaks in zonally-averaged export in the North Pacific, the tropical Pacific and Indian and in the Southern Ocean and low export fluxes in the subtropics and the tropical Atlantic (Fig. 11).

Flux measurements as additional constraint
This is in agreement with the results of Jin et al. (2006) and Lee (2001), with the exception that Lee (2001) Sarmiento and Gruber (2006), and this study. The relatively low North Atlantic export suggested by these latter studies appears to be in conflict with the occurrence of coccolithophorid 640 blooms in this region (Brown and Yoder, 1994). On the other hand, the low TA * inventory in the Atlantic (Fig. 1) argues for a limited export of CaCO 3 in this basin. despite this large riverine source. They conclude that these large riverine inputs must therefore be more than balanced by strong net CaCO 3 formation. For comparison, our estimates of CaCO 3 650 export between 35 • S and 30 • N in the Atlantic are 0.06 (0.018-0.1) Gt-C yr −1 . We note, however, that a direct comparison between these two studies remains difficult. The river input does not strictly have to be compensated by the export flux but rather by the burial flux. Also, the river input gets mixed and is subducted and transported southward by North Atlantic Deep Water and interpretation of concentrations without explicit consideration of transport and mixing is always difficult.

655
Regional patterns of CaCO 3 export (Fig. 10) vary among the different studies. In our approach CaCO 3 export fluxes are scaled by a rain ratio to simulated export of particulate organic carbon within each of the eight considered regions. Thus, the total export for the three Southern Ocean sectors, the Indian Ocean, the tropical and northern Pacific and the tropical and northern Atlantic is constrained by the TA * data, but not the pattern within each of these regions.

CaCO 3 dissolution in the upper ocean
The CaCO 3 leaving the surface ocean dissolves within the water column, at the sea floor or gets buried. Berelson et al. (2007) get to global dissolution rates of ∼ 1 Gt-C yr −1 within upper level waters, clearly higher than our TA * -based estimate of 0.35 (0.26-0.46) Gt-C yr −1 . On a regional level, only two out of ten regional estimates by Berelson et al. (2007) are within our Latin-Hypercube 665 uncertainty ranges; these are the estimates for the tropical Pacific (taken as 5 • N to 5 • S) and the low-and mid-latitude Atlantic (40 • S to 40 • N) ( Table 5). High dissolution rates in the range of the measured decrease of suspended CaCO 3 particles with depth multiplied with a CaCO 3 particle settling velocity of 80 m day −1 and neglecting any temporal trend in the CaCO 3 particle concentration. These estimates may be affected by uncertainties in the assumed particle settling velocity.
CaCO 3 particles settling velocities are reported by Jansen et al. (2002), to vary greatly (0.15 to 3440 m day −1 ) and to be typically order one m day −1 for coccolitophorides and several 100 m day −1 675 for foraminifera and pteropods. Our dissolution rates would be consistent with the measured depth gradient in suspended biogenic CaCO 3 particles for an average settling velocity of a few m day −1 .
As mentioned above, we link the differences between the estimates of Berelson et al. (2007) and this study to methodological problems (Friis et al., 2006) associated with the TA * -CFC age method that very likely introduce a high bias in the results of Berelson et al. (2007).

680
The TA * -CFC age method relies on deduced, observation-derived TA * concentrations and estimates of water mass age, typically derived from measurements of chlorofluorocarbons (CFCs) and their known atmospheric history. TA * concentrations are plotted against their CFC-age and a line is fitted to this data. The higher the TA * concentration for a given water mass, the more TA * must have been added by dissolution to this particular water parcel according to this method. The slope 685 of the relationship between TA * and age is in this sense the CaCO 3 dissolution rate (mol volume −1 time −1 ).
The method has been criticised for its neglect of explicit transport and mixing processes. In particular, Friis et al. (2006) noted that TA * signals ended up above the saturation horizon in their model run, even though there was explicitly no dissolution allowed to occur there. This is confirmed by our 690 sensitivity simulations with no dissolution above the saturation horizon (Fig. 2). This finding does not depend on the choice of the numerical model as mixing within the ocean must spread the TA * signal within the ocean and establish a surface-to-deep gradient in TA * in the upper ocean even when all CaCO 3 dissolves at great depth. The TA * -CFC age method does not account for such processes and assigns dissolution rates in the waters above the saturation horizon irrespective of whether or not 695 the signal stems from the deep ocean. Our approach to combine TA * data within an ocean transport model and the approach by Jin et al. (2006) or Lee (2001) avoid this shortcoming.
In addition, we notice, that in most regions the upper ocean TA * distribution remains remarkably similar within our ensemble, even for very different dissolution rate profiles (see Fig. 2). And in particular, that, for a given export, k slow dissolution rate profiles, with most of the dissolution at or near 700 the ocean bottom, tend to reach highest values of TA * both in the deep ocean and the thermocline (see also inventories in Table 2). In these cases, most of the TA * source is added to slowly ventilated waters in the deep where it accumulates over time along deep water flow paths. This high TA * signal is being brought to the thermocline and eventually to the surface where TA * is reset to its preformed value of zero. As a result, a larger TA * gradient is established across the thermocline for deep com-705 pared to shallow dissolution. The TA * -CFC age method would therefore have a tendency to assign higher dissolution rates in the wrong cases. This further illustrates the difficulty to uniquely relate upper ocean tracer concentrations to either dissolution or mixing processes.

CaCO 3 dissolution in the deep ocean
Turning to the deep ocean, Berelson et al. (2007)  121 Gt-C yr −1 . We note that CaCO 3 formation by coral reefs and burial in shallow coastal waters (Milliman, 1993;Vecsei and Berger, 2004) is not considered in our coarse resolution model.

Parameterisations of CaCO 3 within Earth System Models
Uncertainties appear too large to objectively determine well-defined parameter ranges for the dissolution rates of calcite and aragonite (judging from both the TA * and the flux data compilation).

725
A good agreement between simulated and observed TA * fields can be achieved irrespective of whether dissolution is assumed to depend on the calcite or aragonite saturation state or whether dissolution rates are assumed to be constant throughout the water column. We recall that the computation of the saturation state by carbonate chemistry routines across all grid cells and for each model time step poses a considerable computational burden. For simplicity and to minimise com-730 putational costs, we therefore recommend to describe CaCO 3 dissolution by a constant dissolution rate in Earth System Models as long as these uncertainties exist. This yields an exponential particle flux profile when assuming constant settling velocities throughout the water column. This approach is used in previous studies (Archer and Maier-Reimer, 1994;Ridgwell et al., 2007). A shortcoming of the application of an exponential particle flux profile for CaCO 3 is that it is not easy to account 735 for the potential influence on dissolution of changes in environmental variables, including a decrease in saturation state as expected under ongoing ocean acidification, or in the quality, form, and size distribution of exported CaCO 3 particles. gest that about 30% of the deposition flux dissolves back into the ocean and 70% gets buried in consolidated sediments. 755 We find that the higher the export fluxes within the constrained, likely ranges, the more likely dissolution above the saturation horizon is needed to distribute TA * skillfully in the model. Different kinds of dissolution schemes (with and without dissolution above saturation) achieve realistic TA * distributions within the export ranges identified. Therefore, background dissolution above the saturation horizon cannot be ruled out from this Latin-Hypercube ensemble evaluated within the 760 Bern3D EMIC. Future progress likely depends on a better observational characterisation of particle concentrations, size distribution, and settling within the water column. Physical transport and mixing reduces concentration gradients such that conceptually different dissolution schemes (e.g. dissolution permitted above saturation or not) cannot be distinguished statistically. This implies that concentrations cannot be used to infer dissolution rates directly. It also implies that dissolution pa-765 rameterisations remain uncertain. For simplicity and to minimise computational costs, we suggest to use saturation-independent parameterisations of CaCO 3 dissolution within Earth System Models.

Appendix: Model evaluation and additional constraints
It is an essential prerequisite for the model to feature the main water masses and mixing time scales of the ocean to realistically simulate biogeochemical tracers such as TA * . In this Appendix, we graph-   Table 2. Results from sensitivity simulations. CaCO3 export and TA * inventories for different physical mixing (diapycnal mixing coefficient, kdia) and CaCO3 dissolution schemes. For these illustrative simulations, calcite and aragonite particles were assigned equal parameter values and 10 % of export is assumed to be in the form of aragonite. Fast: k0 = 10 day −1 , n = 1, fcalc = 0.9, kbg = 0; slow: k0 = 0.16 day −1 , n = 2, fcalc = 0.9, kbg = 0; constant: kbg/v = 1/2900 m −1 .          with a higher CaCO3 export tend to require a higher background dissolution to achieve a good skill. [mmol-C m -2 yr -1 ] Figure 10. Observationally-constrained CaCO3 export fields as estimated by different studies. Global total export is estimated to 0.5 by Sarmiento and Gruber (2006), to 1.14 by Jin et al. (2006), to 1.1 by Lee (2001), and to 0.90 Gt-C yr −1 by this study.  Sarmiento & Gruber 2006Jin et al. 2006Lee 2001 This study (median) Figure 11. Zonally-integrated export fluxes of the constrained Bern3D-model ensemble (median and 68 % confidence interval) compared to estimates by Gruber (2006), Jin et al. (2006) and Lee (2001).      Figure A6. As Fig. 9 but coloured according to model skill with respect to the sediment trap database evaluated by basin. This target variable constrains export and dissolution to wider ranges (yellow to red colours) as compared to the TA * target (light green shading). Generally, high skill scores with respect to regional, sedimentcorrected TA * are also associated with high skill scores with respect to regional fluxes. Only few models are in good agreement with the few flux measurements in the Indian Ocean.