2.1 Data

2.1.1 Sediment data and calculations

Sedimentary alkalinity generation due to S burial in the Baltic Proper (sub-basins 7–9 in Fig. 1) was estimated using published data of S contents at F80, a 191 m deep site in the Fårö Deep of the Gotland Sea (Fig. 1; Lenz et al., 2015b). Concentrations of S in micromoles per gram (µmol g⁻¹) were first converted to units of micromoles per cubic centimeter (µmol cm⁻³) using measured porosities and a sediment density of 2.65 g cm⁻³, which is typical for such sediments, and were subsequently depth-integrated between 0 and 25 cm in sediment depth. Following the age model presented by Lenz et al. (2015b), this depth interval represents the burial since 1970, allowing the calculation of an annually averaged rate of S burial (mmol m⁻² yr⁻¹). Using a 1:2 ratio between S burial and TA generation, the latter was calculated and subsequently extrapolated to the basin scale using the total muddy sediment area for the Baltic Proper (Table 1; Al-Hamdani and Reker, 2007).

Table 1Total sediment areas and muddy sediment areas (1000 km²). The muddy sediment areas are based on Al-Hamdani and Reker (2007). Sub-basins according to Fig. 1.

Download Print Version | Download XLSX

Although total S concentrations do not indicate in which form S is buried, this does not matter for the associated TA generation. The conversion from ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ to reduced sulfur produces 2 mol of TA (in the form of ${\mathrm{HCO}}_{\mathrm{3}}^{-}$) per mole of S (Reaction R1, CH₂O represents simplified organic matter). This is irrespective of whether it ultimately ends up in the form of Fe monosulfides (FeS), pyrite (FeS₂), or elemental sulfur (S⁰) when being converted to a solid form. Reductive dissolution of Fe oxides also produces 2 mol of TA (as ${\mathrm{HCO}}_{\mathrm{3}}^{-}$) per mole of dissolved iron (Fe²⁺) formed, but this is compensated for when Fe²⁺ subsequently reacts with dihydrogen sulfide (H₂S) during FeS formation, thereby releasing protons (Reaction R2). Therefore, there is no net TA generation associated with the formation of Fe²⁺ and its subsequent burial as Fe sulfide minerals (Hu and Cai, 2011a).

$$\begin{array}{ll}\text{(R1)}& {\displaystyle }& {\displaystyle }\mathrm{2}{\mathrm{CH}}_{\mathrm{2}}\mathrm{O}+{\mathrm{SO}}_{\mathrm{4}}^{=}⟶\mathrm{2}{\mathrm{HCO}}_{\mathrm{3}}^{-}+{\mathrm{H}}_{\mathrm{2}}\mathrm{S}{\displaystyle }& {\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{\mathrm{CH}}_{\mathrm{2}}\mathrm{O}+\mathrm{Fe}(\mathrm{OH}{)}_{\mathrm{3}}+{\mathrm{H}}_{\mathrm{2}}\mathrm{S}+{\displaystyle \frac{\mathrm{7}}{\mathrm{4}}}{\mathrm{CO}}_{\mathrm{2}}\\ \text{(R2)}& {\displaystyle }& {\displaystyle }⟶\mathrm{FeS}+{\displaystyle \frac{\mathrm{3}}{\mathrm{4}}}{\mathrm{H}}_{\mathrm{2}}\mathrm{O}+\mathrm{2}{\mathrm{HCO}}_{\mathrm{3}}^{-}+\mathrm{2}{\mathrm{H}}^{+}\end{array}$$

2.2 Model calculations

2.2.1 Sediment reactive-transport model (RTM)

A one-dimensional RTM (Reed et al., 2016) was used to calculate benthic TA generation and release at F80, with a minor modification: the redox reaction equations as presented in Table S8 by Reed et al. (2016) were updated to include the total concentrations of the ammonium and sulfide acid–base systems, instead of the corresponding acid–base species. The model calculates acid–base speciation using the direct substitution approach (Hofmann et al., 2008), in which pH and total quantities (dissolved inorganic carbon (DIC), ΣH₂S, total ammonium ($\mathrm{\Sigma }{\mathrm{NH}}_{\mathrm{4}}^{+}=\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\right]+\left[{\mathrm{NH}}_{\mathrm{3}}\right]$), etc.) are used as state variables, meaning that TA is calculated as an output variable. Effluxes of TA from the sediment were subsequently calculated from the gradient in the diffusive boundary layer (Boudreau, 1997). The model includes the carbonate, sulfide, ammonium, and phosphate acid–base systems, such that TA is defined as

$$\begin{array}{ll}{\displaystyle }\mathrm{TA}& {\displaystyle }=\left[{\mathrm{HCO}}_{\mathrm{3}}^{-}\right]+\mathrm{2}\left[{\mathrm{CO}}_{\mathrm{3}}^{=}\right]+\left[{\mathrm{HPO}}_{\mathrm{4}}^{=}\right]+\mathrm{2}\left[{\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}\right]\\ \text{(1)}& {\displaystyle }& {\displaystyle }+\left[{\mathrm{NH}}_{\mathrm{3}}\right]+\left[{\mathrm{HS}}^{-}\right]-\left[{\mathrm{H}}^{+}\right].\end{array}$$

The RTM thus ignores contributions of the fluoride, borate, and silicate acid–base systems, as well as the hydroxide ion and phosphoric acid, which make up part of the classical definition of TA (Dickson, 1981) but were expected to have a low contribution to TA in this setting. Moreover, the model neglects organic alkalinity, which can substantially contribute to Baltic Sea pore water TA (Łukawska-Matuszewska, 2016; Łukawska-Matuszewska et al., 2018), but which is challenging to calculate due to the variety of acid–base groups associated with organic matter. Further details on the governing equations, redox and equilibrium reactions, reaction parameters, and boundary conditions are given in Tables S1–S2 and S4 (Supplement) and in Reed et al. (2016).

The RTM has previously been used at this location to assess the impact of shelf-to-basin Fe shuttling on the formation and stability of the Fe(II) phosphate mineral vivianite (Fe₃(PO₄)₂•8H₂O; Reed et al., 2016). To this end, it was calibrated against a wide selection of pore water and solid phase data presented in Jilbert and Slomp (2013), Lenz et al. (2015b), and Reed et al. (2016). We confirm this calibration for the carbonate system with additional DIC and TA pore water data from a multicore recovered at F80 during a research cruise with R/V Pelagia in June 2016. Core handling and pore water analyses have been performed following Egger et al. (2016). We used the previous calibration and perform sensitivity analyses to identify the key mechanisms responsible for benthic TA generation and release. To represent the variety of bottom water conditions at site F80 since the 1970s, four time intervals were recognized (Fig. 2): 1970–1973 (baseline; I), 1973–1978 (start of change in Fe loading; II), 1978–1981 (eutrophication but pre-euxinia; III), and 1981–2009 (eutrophication and euxinia, IV).

Figure 2Variability in bottom water redox conditions, organic carbon, and Fe inputs between 1970 and 2009, used to force the reactive-transport model at site F80. Numbers indicate the four different time intervals recognized in this study: (I) baseline (1970–1973), (II) start of change in Fe loading (1973–1978), (III) eutrophication but pre-euxinia (1978–1981), and (IV) eutrophication and euxinia (1981–2009). Figure modified from Reed et al. (2016).

Download

3.1 Sediment and RTM calculations

Both the model- and observation-based estimates indicate that between 1970 and 2009, on average 291–295 mmol S m⁻² has annually been buried, leading to a TA generation of 582–590 mmol m⁻² yr⁻¹. This corresponds to a TA flux of ∼43.2–43.8 Gmol yr⁻¹ from muddy sediments (Table 2). The model further suggests that virtually all of the S solids are present in the form of FeS₂. Of the 291 mmol m⁻² yr⁻¹ of S being buried, only 56 % was formed in situ, whereas the remaining 44 % was deposited as a result of the shuttling of Fe in the form of FeS₂ to the deep basin (Lenz et al., 2015a). However, as BALTSEM does not resolve FeS₂ formation in either the water column or sediment, both need to be included when estimating the unresolved TA source due to ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction and S burial.

Table 2Estimated S burial (mmol m⁻² yr⁻¹) and associated TA generation (Gmol yr⁻¹) for site F80 (lat 58.0000, long 19.8968) between 1970 and 2009. Both S contents and dating were taken from Lenz et al. (2015b). The basin-scale calculation was based on the total muddy sediment area for the Baltic Proper in BALTSEM of 74 300 km² (Table 1). Numbers in parentheses represent S burial due to in situ S formation only, thus excluding S burial due to FeS₂ deposition.

Download Print Version | Download XLSX

In line with other studies (e.g., Jørgensen, 1977), the vast majority of reduced S produced through ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction was reoxidized in either the sediment or the overlying water. On average, only 10.2 % was buried, but there was strong temporal variability in this percentage (Table 3). The fraction of S solids being buried was highest under eutrophic but non-euxinic conditions (i.e., period III; 40.7 %). Since 1981 (period IV), inputs of Fe oxides have decreased, leading to a higher efflux of ∑H₂S and thus less S burial, despite higher ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction rates (SRRs). Our results indicate that even under non-euxinic conditions, pyrite formation was limited by the availability of highly reactive Fe, as is the case in most marine systems (Berner, 1984; Raiswell and Canfield, 2012). This limitation was confirmed by the difference between potential and simulated S formation rates (Table 3), for which the former indicates the amount of solid S that could have formed based on the other modeled sources and sinks of ∑H₂S. It is thus indicative of the amount of S mineral formation under unlimited Fe²⁺ supply.

Table 3Estimated depth-integrated ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction (SRR), S reoxidation (S-OX), and S⁰ disproportionation rates (S⁰⁻dispr), as well as ∑H₂S efflux (all in mmol S m⁻² yr⁻¹), as derived from the one-dimensional reactive-transport model (Reed et al., 2016) for the periods 1970–1973 (baseline; I), 1973–1978 (start change in Fe loading; II), 1978–1981 (eutrophication but pre-euxinia; III), and 1981–2009 (eutrophication and euxinia; IV). The difference between SRR (production of reduced S) and the sum of S-OX, S⁰⁻dispr, and ∑H₂S efflux (removal of reduced S) is assumed to be representative for the maximum potential S formation. The simulated S formation is derived from the mass balance of the RTM (see Table S6 for details).

Download Print Version | Download XLSX

Figure 3Pore water profiles of selected simulated (solid lines) and observed (dots) variables at site F80 (a–h); simulated rates of some major processes impacting TA dynamics at the end of the four major time intervals (i–p; all rates in mmol TA dm⁻³ yr⁻¹). The blue area in (m) indicates the sulfate–methane transition zone (SMTZ) in 2009; in the other years, the SMTZ was located around the depth of the modeled interval (32 cm). Additional pore water and solid phase profiles were published in Reed et al. (2016). Previously unpublished measurements can be found in Table S5 (Supplement). Dates in (a)–(h) indicate sampling dates; LC: long core.

Download

In the period 1970–2009, S burial could on average only explain 328 mmol m⁻² yr⁻¹ of TA generation (Table 3; using a 1:2 ratio between S burial and net TA generation). This was 9.2 % of the internally generated TA (3948 mmol m⁻² yr⁻¹), but again clear temporal variations were observed (Table 4). Under the baseline conditions (period I), when little S was buried, it only made up 3.7 % of the total TA generation. This percentage increased to 12.4 % between 1973 and 1978 (period II), when more Fe was available, and peaked at 34.5 % between 1978 and 1981 (period III), when both Fe(OH)₃ and carbon loadings were high. Since 1981 (period IV), the decrease in Fe(OH)₃ loading has further limited S burial, leading to a contribution of only 6.7 % of the internally generated TA.

Pore water profiles of DIC and TA (Fig. 3f, g) indicate that the model is generally well calibrated for the carbonate system. Both DIC and TA concentrations were lower in 2016 compared to the earlier measured data against which the model was calibrated (Reed et al., 2016). The profiles furthermore show that the model overestimates the buildup of ∑H₂S (Fig. 3c). This can be explained by loss of ∑H₂S during sampling, which is a common problem for anoxic sediments, but also by the chosen lower boundary condition of the model. Whereas the model assumes no gradient with the underlying sediment, the data for 2009 (Jilbert and Slomp, 2013; Lenz et al., 2015b) show a declining trend with depth below 32 cm, suggesting a downward diffusive flux of ∼144 mmol ∑H₂S m⁻² yr⁻¹ in 2009. When fixed at depth as Fe sulfides, this flux would lead to an additional TA generation of ∼288 mmol m⁻² yr⁻¹ for 2009. However, given the lack of major S accumulation in the sediment below 32 cm (Jilbert and Slomp, 2013), we do not believe this downward flux contributed greatly to net TA generation over the full period of investigation.

Table 4Estimated TA generation from the reactive-transport model integrated over the whole sediment column for the periods 1970–1973 (baseline), 1973–1978 (start change in Fe loading), 1978–1981 (eutrophication but pre-euxinia), and 1981–2009 (eutrophication and euxinia) for the dominant processes, as well as TA generation and efflux (all in mmol m⁻² yr⁻¹).

Download Print Version | Download XLSX

Rate profiles of the most important processes contributing to TA (Fig. 3i–p). show that especially between 1978 and 1981 (period 3), when OM and Fe inputs were high but bottom waters were still oxic, intense cycling of Fe occurred in the sediments, associated with high TA production and consumption. Dissolved Fe produced from reductive dissolution of amorphous Fe oxides during OM degradation either diffused upward, where it reoxidized in the oxic sediments (Fig. 3n), or downward to react with ∑H₂S (Fig. 3o). Well-crystallized Fe oxides, assumed to be inaccessible for OM degradation, reacted with ∑H₂S over a wide range, thereby producing additional Fe²⁺ (Fig. 3p).

On a system scale, this cycling of Fe does not lead to net TA generation (Hu and Cai, 2011a), but it may impact the efflux of alkalinity from the sediment that is calculated by the RTM. This flux cannot directly be used to assess the long-term net TA generation that we are interested in, as it is the product of a variety of reversible and irreversible TA-generating reactions, such as the intense Fe cycling discussed above. Moreover, its constituents (e.g., HS⁻) may become reoxidized in the water column. However, the magnitude and temporal variability in the efflux compared to those of S burial and total TA generation may provide information on its driving processes. Note that the difference between total TA generation and efflux (Table 4; Fig. 4) reflects the buildup of TA in the sediment, as well as loss of TA at depth through burial.

A comparison of their temporal variabilities shows that the benthic TA efflux only partly followed the pattern in S burial (Fig. 4; Table 3). Since 1973, when the efflux was 1901 mmol m⁻² yr⁻¹, it decreased to 1561 mmol m⁻² yr⁻¹ in 1978, followed by a sharp increase to 3261 mmol m⁻² yr⁻¹ in 1982 and a more gradual increase to 4823 mmol m⁻² yr⁻¹ in 2009. Generation of TA throughout the entire sediment column contributed to the calculated efflux (Table 4), to a major extent due to high rates of ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction at depth (Fig. 3k, m). The methane diffusing upward from deeper sediment layers played a key role here, being responsible for on average ∼95 % of the CH₄-driven ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction and ∼43 % of the total ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction. The temporal change in spatial pattern of the CH₄-driven ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$e reduction (Fig. 3m) is a direct result of more ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ being consumed by organic matter degradation since the start of bottom water euxinia (Fig. 3k, Table 4). This is confirmed by the upward-shifting sulfate–methane transition zone (SMTZ) since the onset of bottom water euxinia (Fig. 3d; see also Reed et al., 2016), whose position matches the highest reaction rates of CH₄-driven ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction (Fig. 3m). The minor peak in CH₄-driven ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction around 5 cm in depth in 1973 and 1978, and the slightly more pronounced peak at ∼3 cm in depth in 1981, are in contrast driven by in situ produced methane due to methanogenesis (Fig. 3l).

Figure 4Calculated TA efflux and generation at site F80 integrated over the whole sediment column (all in mmol m⁻² yr⁻¹) due to various biogeochemical processes implemented in the RTM.

Download

Interestingly, the decrease in efflux between 1973 and 1978 (period II), which resulted from changes in the Fe loading, was not mimicked in either the total TA generation or the amount of S burial. Rather, it reflected the pattern of the change in TA generation through secondary reactions (Fig. 4). The most important secondary reaction contributing negatively to TA between 1973 and 1978 was the reoxidation of Fe²⁺, the rate of which more than doubled during this period (Table 4) and which, in contrast to the other dominant reactions, was restricted to the upper centimeter of the sediment column (Fig. 3n). This indicates that it was the driving force of the lower TA efflux during this period. Although reoxidation of Fe²⁺ consumed even more TA in period III (1978–1981), this was more than compensated for by the concurrent enhanced TA generation due to OM degradation, especially coupled to Fe-oxide reduction, even though that occurred deeper in the sediment (Fig. 3j; Table 4).

Table 5Average resolved and unresolved TA sources minus sinks (SMS) and river loads with standard deviations (Gmol yr⁻¹) in 1970–2014 according to the BALTSEM calculations in this study. Sub-basins according to Fig. 1.

Download Print Version | Download XLSX

Figure 5Annual mean observed (dots) and modeled (lines) normalized surface water TA (TA_N, µmol kg⁻¹) and salinity in five sub-basins. Full and dashed lines represent the scenarios in which the unresolved TA sources are added as land loads or sediment release, respectively. Model data from sub-basins 3 (SK), 9 (GS), 10 (BS), 11 (BB), and 13 (GF) are compared to observed data at the Anholt East, BY15, US5B, F9/A13, and LL7 stations, respectively.

Download

In summary, this discrepancy between sedimentary TA generation due to S burial and modeled effluxes of TA highlights that both represent processes acting at various spatial and temporal scales. Long-term TA generation should be interpreted as the net TA generation, i.e., the TA change occurring after all reoxidation reactions took place, in the coupled water column–sediment system. In contrast, calculated efflux of TA, as well as TA generation through various processes at a specific moment in time within different zones in the sediment, is highly variable and is directly impacted by local coupled dynamics of S, Fe, and CH₄ (see also Table 4, Fig. 3). While the sedimentary processes are highly relevant to understanding the major factors driving short-term TA dynamics, ultimately it is the burial of S that represents a TA source relevant in the long term (see also Sect. 2.2.3).

3.2 BALTSEM calculations

The recalibrated unresolved TA sources as well as the resolved pelagic and benthic TA sources minus sinks in the different sub-basins as calculated with BALTSEM are indicated in Table 5. In total, the recalibrated unresolved source amounts to 260±0 Gmol yr⁻¹ (constant over time), while the total resolved pelagic and benthic sources minus sinks amount to a net source of 120±47 Gmol yr⁻¹ over the 1970–2014 period. For comparison, the riverine TA load amounts to 470±62 Gmol yr⁻¹. In Figs. 5–6, simulated and observed surface and deep water TA normalized to mean salinity at the corresponding station and water depth (TA_N) and salinity are shown. The normalized TA is used in order to avoid uncertainties related to discrepancies between simulated and observed salinity. Full lines in Figs. 5–6 represent the scenario in which the unresolved sources were added as land loads, whereas dashed lines represent the case in which these sources were instead modeled as sediment effluxes.

Figure 6Annual mean observed (dots) and modeled (lines) normalized deep water TA (TA_N, µmol kg⁻¹) and salinity in five sub-basins. Full and dashed lines represent the scenarios in which the unresolved TA sources are added as land loads or sediment release, respectively. Model data from sub-basins 3 (SK), 9 (GS), 10 (BS), 11 (BB), and 13 (GF) are compared to observed data at the Anholt East, BY15, US5B, F9/A13, and LL7 stations, respectively.

Download

The temporal development of resolved and unresolved TA sources minus sinks throughout the model simulation is shown in Fig. 7. In this simulation, the unresolved sources in the different sub-basins were assumed to remain constant throughout the model run (Fig. 7), while the resolved sources and sinks vary depending on primary productivity, oxygen conditions, denitrification rates, and other biogeochemical processes included in the model. Despite the constant unresolved sources, simulated TA_N concentrations generally reproduce observed values. Exceptions are the simulated TA_N concentrations in the Kattegat and the Gotland Sea in the 1980s, where actual concentrations are overestimated, and TA_N concentrations in the Bothnian Sea and Bay that are underestimated in the last 10-year period (Figs. 5–6).

There is an overall long-term increase in the resolved net TA generation in sediments and water column combined (Fig. 7), reflecting the ongoing eutrophication and overall deteriorating oxygen conditions of the Baltic Sea. The resolved net pelagic TA source increases in the period 1970–2000 in response to an increased primary production and then levels out and slightly declines in the last decade. The increased resolved benthic source in the last decade, however, is a response to deteriorating oxygen conditions resulting in increased TA generation through denitrification and ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction. Sulfate reduction in the BALTSEM model is however not an irreversible source since sulfidic waters can be reoxidized by deep water inflows, thus consuming TA and reversing the source.

4.1 Use of BALTSEM and RTM in the context of this work

Given the detailed presentation of sedimentary processes and effluxes in Sect. 3.1, one may wonder why only S burial is used in the coupling to BALTSEM. After all, the RTM calculations include many processes other than S burial. However, to study the impact of sedimentary TA generation on the long-term TA development in the Baltic Sea, we need to take into account only those processes that are relevant to accomplish this task.

BALTSEM includes many biogeochemical processes that produce and consume TA both reversibly and irreversibly on short timescales and in many boxes within each sub-basin of the Baltic Sea. These processes are described in Sect. 2.2.2 and are further listed in detail in Table S3. BALTSEM furthermore accounts for land loads, atmospheric depositions, and TA exchange between sub-basins and between the Baltic Sea and the North Sea. The result of the model simulations, i.e., the long-term development of TA in various sub-basins, is what we compare to observations in the water column (Figs. 5–6). Similarly, the RTM calculates net TA generation due to various reversible and irreversible processes (described in detail in Tables S1–S2). If we dynamically coupled the RTM to BALTSEM, we would have to consider all these processes and link all species between both models. Given the unfeasibility of this, as discussed in Sect. 2.2.3, we couple both models by using the output of the RTM to further constrain BALTSEM. Specifically, we explain part of the source of BALTSEM that is unresolved but necessary to describe the long-term TA development in the Baltic Sea.

Figure 7Annual mean TA sources minus sinks (SMS) (Gmol yr⁻¹) in the entire Baltic Sea according to BALTSEM calculations: Resolved benthic sources minus sinks (SMS) (green line), resolved pelagic SMS (purple line), total resolved SMS (yellow line), unresolved SMS (red line), and total SMS (blue line).

Download

This means that in this context we only need to consider the processes from the RTM that are (a) irreversible on the timescale of interest (i.e., decades) and (b) not included in BALTSEM. Burial of Fe sulfides (Hu and Cai, 2011a) is the only major process that falls in this category. Denitrification using an external ${\mathrm{NO}}_{\mathrm{3}}^{-}$ source, the other main pathway for net TA generation (Hu and Cai, 2011b), is already included in BALTSEM. Many other sedimentary processes produce or consume TA (Tables S1–S2, Supplement; Soetaert et al., 2007), but they are not irreversible on the relevant timescale. Their dynamics are, however, highly interesting to discuss as they help determine what limits net sedimentary TA generation and which processes mainly drive the effluxes of TA and other constituents to the water column. Note that this irreversibility is also a reason why we do not use these effluxes as input to BALTSEM. In addition, they are already partly included in BALTSEM, e.g., in the case of ∑H₂S produced from ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction.

The RTM fluxes are upscaled under the assumption that the fluxes computed for the F80 site are representative for the muddy sediment area of the Baltic Proper. This assumption is associated with uncertainties because of spatial differences in the sediment geochemistry of muddy Baltic Proper sediments as illustrated by the pore water and Fe–S chemistry for four other sites as published by Lenz et al. (2015). The solid phase profiles for these sites show temporal trends over the past decades similar to those of F80. Furthermore, the pore water profiles show that site F80 has a relatively high rate of organic matter deposition and alkalinity regeneration when compared to most of the other sites. This implies that, with our extrapolation, the role of the sediment could be slightly overestimated. Thus, the large-scale fluxes we obtain by extrapolating fluxes from one specific site are to be regarded as a maximum estimate of the contribution of S burial to the overall TA budget of up to 26 %.

Although a full coupling between the two models is not a realistic goal at the moment, the development of sediment processes in BALTSEM is decidedly a highly desirable future goal. In particular, the inclusion of sedimentary Fe–S dynamics and related phosphorus (P) cycling would serve to improve our understanding of both TA and P dynamics on a system scale. The present study can be seen as an intermediate step towards a more detailed (if not complete) model description of sediment processes in the Baltic Sea. In fact, the relatively large influence of sedimentary processes on TA dynamics that we demonstrate in this study also serves as a motivation to pursue this goal.

4.2 Sulfur burial and TA generation in the Baltic Proper

While mineralization in the sediments occurs everywhere where there is labile organic matter, permanent burial of organic matter as well as other solids such as Fe sulfides should predominantly occur in muddy sediments. Consequently, the part of the unresolved TA source that is a result of S burial should be released from muddy sediments rather than from the entire sediment surface area of the Baltic Sea. Our RTM calculations in combination with observations from site F80 in the Baltic Proper (sub-basins 7–9 in Fig. 1) provide the first estimate of TA generation resulting from S burial in Baltic Sea sediments (582–590 mmol m⁻² yr⁻¹). Assuming that the calculated TA generation resulting from S burial is representative only for the muddy sediment area in the Baltic Proper (74 300 km²; Table 1), the total annual flux in this area is ∼44 Gmol yr⁻¹. The calibrated unresolved TA source in the Baltic Proper amounts to 170±0 Gmol yr⁻¹ according to the BALTSEM model (Table 5).

In the two different scenarios in which the unresolved source is added as land loads (full lines in Figs. 5–6) or sediment release (dashed lines in Figs. 5–6), the simulated surface water TA concentrations are very similar (Fig. 5). Deep water concentrations, however, differ significantly in the Gotland Sea and the Gulf of Finland but not in the other sub-basins (Fig. 6). The reason behind the rather similar results for these two different scenarios is that land loads supplied to the different basins are rapidly distributed in the well-mixed surface layer, and the well-mixed surface layer constitutes a large majority of the water volume. In the deeper and more isolated parts of the system, TA concentrations are lower in the “land loads” case compared to the “sediment release” case.

In the sediment release case, the unresolved TA source in the Baltic Proper (170±0 Gmol yr⁻¹) corresponds to a flux of 730 mmol m⁻² yr⁻¹ if the source is distributed evenly over the entire sediment surface (228 000 km²). However, if the unresolved source is instead constrained only to muddy sediments, the flux would amount to 2236 mmol m⁻² yr⁻¹, which is far above the long-term mean flux due to S burial as obtained by RTM calculations. Even during peak pyrite formation periods, S burial only resulted in a source of 1078 mmol TA m⁻² yr⁻¹ (Table 3). Furthermore, in a BALTSEM experiment in which the unresolved sources were released only from muddy sediments (but at higher rates corresponding to the smaller surface areas), the deep water TA concentrations in particular in the Baltic Proper were overestimated while the surface water TA concentrations were underestimated (not shown). Based on the RTM calculations, TA generation coupled to S burial could thus account for 26 % of the unresolved source, at least in this sub-area of the system. The remaining unresolved TA source of ∼74 % could possibly be explained by underestimated river loads or submarine groundwater discharge of TA (e.g., Szymczycha et al., 2014). We have no data to quantify these fluxes, however.

4.3 Reversible versus irreversible sedimentary processes generating TA

As demonstrated by both simulated and observed sediment profiles at F80, a transition from hypoxic to euxinic conditions around 1980 resulted in a strong increase in both solid phase S and Fe burial (Reed et al., 2016). Furthermore, the molar S-to-Fe ratio of ∼2 suggests formation and burial of mostly FeS₂. Both FeS₂ and FeS can be formed from reactive Fe²⁺ and sulfide produced during Fe(OH)₃ and ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction, respectively. These redox reactions ultimately result in net TA generation (Tables S1–S2, Supplement). Another possible pathway is that methane (formed by methanogenesis) is oxidized anaerobically by reduction of either ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ or Fe(OH)₃ (Slomp et al., 2013; Egger et al., 2015b), and Fe and sulfide can then be sequestered in the form of FeS₂. Results from the RTM indicate that CH₄ and organic matter are both important electron donors at F80. CH₄ oxidation contributes to, on average, 43.8 % of total SRR and occurs at greater depth than ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction through organic matter degradation (Fig. 3k, m). Iron-mediated anaerobic oxidation of methane is not included in the set of reactions of this RTM. Previous work has indicated that this process mainly occurs in organic-poor sediments depleted in ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ (Riedinger et al., 2014; Egger et al., 2017). These conditions are not met at F80, rendering an important role for this process unlikely.

Apart from such eutrophication-induced changes in the coupled Fe–S cycling, increasingly euxinic conditions also influence manganese (Mn) sequestration in sediments. Dissolved Mn²⁺ can be sequestered in the form of Mn carbonates (MnCO₃). If this occurs, the TA generation associated with the reduction of manganese dioxide (MnO₂; 2 mol of TA per mole of MnO₂; Tables S1–S2, Supplement) is completely compensated for by the TA sink associated with carbonate removal. However, under euxinic conditions, manganese sulfide (MnS) can be formed if the sulfide availability exceeds the Fe²⁺ availability (Lenz et al., 2015b). Indeed, long-term sediment records indicate a relation between euxinic periods in Baltic Sea deep waters and burial of Mn sulfides in the forms of both rambergite and alabandite (Lepland and Stevens, 1998). As opposed to burial of MnCO₃, burial of MnS results in a net TA generation comparable to that of FeS₂ burial. In the RTM, we did not investigate the possible impact of MnS formation on TA generation as the sediment record at F80 does not show substantial Mn enrichments in the surface, despite higher sulfide than Fe²⁺ availability (Lenz et al., 2015b).

Ammonium, ΣH₂S, Fe²⁺, and Mn²⁺ are rapidly oxidized if oxygen is supplied to anoxic waters. The result is a TA sink that compensates for the TA generation by anaerobic mineralization. Precipitates such as FeS and FeS₂ can also be oxidized but this is generally a slower process, especially in the case of FeS₂ (Millero et al., 1987; Wang and Van Cappellen, 1996). Moreover, these S minerals are embedded in organic-rich, ΣH₂S-producing sediments. This reduces the impact of possible reoxygenation of the sediment for extended periods of time, implying that the TA source that results from S burial is stable. Sediment cores indicate the presence of Fe sulfides – in particular FeS₂ – in the top 3 m of Gotland Sea deep water sediments (Boesen and Postma, 1988) as well as in the top 10 m of deep Bornholm Basin sediments and the top 27 m of Landsort Deep sediments (Egger et al., 2017), all corresponding to roughly 8000 years. In our model results, FeS₂ is the dominant form of S in the sediment. Our work also shows that reoxidation of reduced S never exceeds 0.5 % between 1970 and 2009, irrespective of whether S solids or the total reduced S pool (i.e., including ∑H₂S) is investigated (Table S5, Supplement).

Vivianite formation is another process that generates TA in a net sense. The presence of vivianite in sediment cores (Egger et al., 2015a) indicates that this mineral can be stable upon burial. However, vivianite dissolves in the presence of ΣH₂S (Dijkstra et al., 2018), excluding its burial to be a long-term TA source at F80. It could however be a stable TA source at locations where ΣH₂S rather than Fe availability limits pyrite formation, such as the Bothnian Sea (Egger et al., 2015a).

Calcifying organisms that build calcium carbonate (CaCO₃) shells have a large influence on the carbonate system in many marine areas, as illustrated by high TA fluxes related to CaCO₃ formation and dissolution in the North Sea (Brenner et al., 2016). CaCO₃ formation results in a TA drawdown in the productive layer and a TA source where the shells are dissolved. Burial of CaCO₃ shells is a net TA sink on a system scale. In the Baltic Sea, however, planktonic calcifiers are largely absent – likely because of low saturation values of calcite and aragonite in winter (Tyrrell et al., 2008). In the RTM, CaCO₃ dissolution and precipitation are included, based on observed sedimentary CaCO₃ contents of ∼200 µmol g⁻¹ (Fig. 3e; see Reed et al., 2016, for further details). The prescribed input of 86 mmol m⁻² yr⁻¹ in combination with prevailing conditions in the sediment led to a net TA loss due to CaCO₃ dissolution, which is on average only −0.01 mmol m⁻² yr⁻¹ (data not shown). For this reason, we have not included the effects of CaCO₃ precipitation and dissolution in our analysis.

4.4 Long-term development of TA in the Baltic Sea

Several studies indicate essentially linear TA–salinity relations in Baltic Sea surface water (Ohlson and Anderson, 1990; Thomas and Schneider, 1999; Perttilä et al., 2006; Beldowski et al., 2010). Long-term TA increases that are not connected to salinity are, however, apparent from observed TA–salinity relations (Fig. 8; Tables S7–S8, Supplement). Furthermore, Müller et al. (2016) found generally increasing TA concentrations decoupled from salinity in the Baltic Sea over the past 2 decades.

Figure 8TA concentrations normalized to salinity = 7 (TA_N), using the TA–salinity relations for the Baltic Proper – Kattegat (full circles), and Baltic Proper – Gulf of Bothnia (open circles) (see Tables S7–S8, Supplement). The Ruppin (1909) value was based on measurements in 1906–1907 (see Dyrssen, 1993), while the Buch (1945) value was based on measurements in 1927–1935 (see Buch, 1945). The Buch (1945) and Perttilä et al. (2006) values were converted to micromoles per kilogram from micromoles per liter.

Download

The resolved pelagic and benthic TA sources minus sinks in the BALTSEM calculations (Fig. 7) on average increase by approximately 3 Gmol yr⁻¹ in the 1970–2014 period. Furthermore, the RTM calculations (Table 3) indicate that S burial can increase by a factor of 4 after a transition from oxic to anoxic/euxinic conditions, and even by an order of magnitude during this transition if both Fe-oxide and organic matter loadings are enhanced. Thus, the fraction of the unresolved source that is a result of S burial should be quite variable depending on mineralization rates and oxygen conditions in different areas of the Baltic Sea as well as during different periods in time.

Anaerobic mineralization occurs in sediments even if the overlying water is oxic, and for that reason TA release coupled to S burial does not exclusively occur from sediments covered by sulfidic waters. In fact, our results indicate the highest S formation rates under eutrophic, but non-euxinic, conditions (Table 3). However, large-scale and long-term changes in TA generation related to changes in S burial are mainly expected to occur in areas experiencing transitions between oxic and anoxic conditions and in addition as a result of changes in Fe loadings (Lenz et al., 2015a). Hypoxic and anoxic conditions in the Baltic Sea water column – as well as rapid transitions between oxic and anoxic conditions – occur primarily in the deep basins of the Baltic Proper, although episodes of oxygen depletion can also occur in the deep water of, in particular, the Gulf of Finland, as well as in many eutrophic coastal fjords and bays.

The long-term TA decrease in the Gotland Sea in the 1980s coincided with a decreasing salinity (Figs. 5 and S1) as well as improved oxygen conditions in large volumes of the deep water (not shown). During this period, stratification was considerably weakened and as a result the halocline depth in the Baltic Proper increased, and inflowing new deep water ventilated primarily the upper deep water. Thus, a much larger water volume than usual was well ventilated (e.g., Stigebrandt and Gustafsson, 2007). It is possible that during this period, S burial and associated TA generation were considerably weakened. A very strong TA increase observed in the early 1990s coincides (more or less) with a strengthened stratification due to saltwater inflows in 1993. A rapid deterioration of oxygen conditions followed because of a suppressed deep water ventilation during periods of strong stratification. This development towards increasingly anoxic/euxinic conditions could potentially cause a large response in terms of TA generation.

It is however likely that the observed TA decline in the 1980s followed by the strong TA increase in the 1990s is exaggerated because of unreliable measurements before 1993. After that, the precision of TA measurements appears to have increased considerably as is evident from the relatively low scatter in TA values after 1993 compared to before 1993 (see Müller et al., 2016, their Fig. 3). In particular, we find in Fig. 6 that even in the Kattegat deep water the observed TA concentrations in the period ∼ 1985–1992 are comparatively low. While the Kattegat surface waters are heavily influenced by outflowing Baltic Proper water, the Kattegat deep water is affected only very marginally. Furthermore, the observed TA values from the Gotland Sea deep water in approximately the same period are considerably lower than our modeled values (Fig. 6).

TA is also influenced by atmospheric deposition on the water surface due to emissions from land and ships (Hassellöv et al., 2013; Hagens et al., 2014). Deposition of S and N oxides represents a TA sink, while $\mathrm{\Sigma }{\mathrm{NH}}_{\mathrm{4}}^{+}$ deposition is a TA source. The net effect is a TA sink; the impact peaked in the 1980s, but has since then diminished due to reduced land emissions (Omstedt et al., 2015). According to our BALTSEM calculations, the TA sink related to acidic depositions has declined from approximately −40 Gmol yr⁻¹ in the 1980s to −10 Gmol yr⁻¹ in the past decade. This reduced TA sink thus contributes to the increasing TA concentrations in the Baltic Sea.

Riverine TA concentrations can increase as a result of enhanced weathering of carbonate and silicate rocks in the catchments. The rate of weathering depends on temperature, precipitation, soil organic matter contents, and deposition of acids (Ohlson and Anderson, 1990; Dyrssen, 1993; Sun et al., 2017). Average TA loads from Swedish rivers in 1985–2012 and Finnish rivers in 2001–2012 amount to 42 and 15 Gmol yr⁻¹, respectively (Fig. S2, Supplement), together corresponding to some 12 % of the total TA load (∼470 Gmol yr⁻¹) to the system. There was a strong long-term increase in the flow-normalized TA loads from Swedish rivers – approximately 21 % over the period 1985–2012. In Finnish rivers, however, there was no increase in the period 2001–2012. Because of a generally poor availability of river data from other countries around the Baltic Sea, we have no clear understanding of the long-term TA development in the great rivers in the southeastern Baltic Sea (where the highest TA concentrations are also generally observed).

Riverine TA concentrations in the BALTSEM model were calculated from observed monthly mean values only in the period 1996–2000. If the long-term increasing trend observed for TA loads in Swedish rivers is also representative for rivers in the southeastern Baltic Sea, this would signify that the model is forced by too low riverine TA concentrations in the last decade but conversely too high concentrations in the first couple of decades. It is plausible that changes in river water properties are responsible for at least part of the overall increasing TA concentrations in the Baltic Sea. This could in particular be the case for the Bothnian Sea and Bay where simulated TA in the last decade is underestimated by the BALTSEM model (Figs. 5–6).

4.5 Simulations of future scenarios

High productivity and deep water oxygen consumption rates favor TA-generating anaerobic mineralization processes. One potential consequence is that a large-scale recovery from eutrophication could reduce the CO₂ buffering capacity of a marine system and thus also reduce the atmospheric CO₂ sink and surface water pH.

In this section we investigate how the simulated TA in BALTSEM responds to two different nutrient load scenarios: (1) the business-as-usual (BAU) scenario with high nutrient loads throughout the 21st century and (2) the Baltic Sea Action Plan (BSAP) scenario with large reductions in N and P loads (Fig. S3, Supplement). We use the ECHAM5 A1B no. 1 scenario for CO₂ emissions and climate change downscaled for the Baltic Sea region (see Omstedt et al., 2012). The A1B emission scenario represents a socioeconomic development producing medium CO₂ emissions for which the atmospheric CO₂ partial pressure (pCO₂) reaches some 700 µatm by the year 2100 (Fig. S4, Supplement). The unresolved TA source is kept constant throughout these simulations. This means that any simulated changes in TA are related to changes in river loads and exchange with the North Sea, as well as changes in TA-producing and TA-consuming biogeochemical processes that are included in BALTSEM (production, mineralization, denitrification, nitrification, ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ reduction, ΣH₂S oxidation, etc.).

Figure 9(a) Simulated annual mean surface water TA and pH in sub-basin 9 (GS) according to the BSAP (black lines) and BAU (red lines) nutrient load scenarios. (b) Differences between the BSAP and BAU scenarios.

Download

According to the BALTSEM simulations, the surface and deep water temperatures in the Gotland Sea will increase by approximately 3 ^∘C over the 21st century, while salinity is reduced by more than 2 (Fig. S5, Supplement). Surface water phosphate concentrations will decline by ∼0.2 µmol kg⁻¹ in the BSAP scenario, resulting in a reduced primary production and increased deep water oxygen concentrations (Fig. S6, Supplement). The reduced productivity and large-scale recovery from anoxic deep water conditions in the BSAP scenario will also have large consequences for TA and in extension CO₂ buffer capacity and pH. Towards the final decades of the simulations, surface water TA in the BAU scenario exceeds that in the BSAP scenario by ∼150 µmol kg⁻¹ (Fig. 9). As a result, the surface water pH is reduced by 0.1 units more in the BSAP than in the BAU scenario.

Table 6Input of Fe oxides and simulated S burial and TA efflux averaged for the period 2011–2050 under a range of environmental conditions as calculated with the reactive-transport model. All values are in mmol m⁻² yr⁻¹.

Download Print Version | Download XLSX

These scenario simulations do not include changes in TA generation resulting from changes in S burial driven by productivity and Fe-oxide availability since these processes are not resolved in BALTSEM. To investigate how the sediment, and more specifically S formation and burial, will respond to changes in Fe and organic carbon loadings, we ran the RTM for an additional 40 years under the present environmental conditions, as well as under a range of changes in these loadings. This sensitivity analysis (Table 6) shows that reverting the productivity regime to pre-1978 conditions decreases the calculated TA efflux by ∼50 %, a direct result of less organic matter degradation, whereas S burial is hardly impacted as it is still limited by Fe availability. Lowering the Fe-oxide loading to pre-1973 values decreases the S burial by an order of magnitude, confirming its limitation by Fe. Strikingly, the TA efflux is only marginally impacted, indicating again the decoupling between short-term flux dynamics and long-term TA generation, as discussed extensively in Sects. 2.2.3 and 3.1. Increasing the Fe-oxide loading to the peak values of 1981 slightly lowers the TA efflux while more than doubling S burial, a direct result of a higher Fe availability. In summary, our sensitivity analysis confirms that the form and rate of Fe input exerts the dominant control on S burial and long-term TA impacts, whereas the rate of organic matter input mainly drives the short-term variability in TA effluxes. It also highlights that sedimentary TA generation due to S burial and modeled effluxes of TA should be regarded as occurring on various temporal scales.

It is a simple exercise to examine the sensitivity of pH to further changes in TA. Using the CO2SYS software (van Heuven et al., 2011; http://cdiac.ess-dive.lbl.gov/ftp/co2sys/CO2SYS_calc_MATLAB_v1.1/, last access: 14 November 2016), and assuming that the surface water pCO₂ is in equilibrium with the atmosphere, a surface water pCO₂ of 700 µatm and TA of 1425 µmol kg⁻¹ (as at the end of the BSAP scenario) results in a surface water pH of 7.80, assuming a salinity of 5.2 and temperature of 11 ^∘C (see Fig. S5). For example, decreasing the surface water TA to 1325 or 1225 µmol kg⁻¹ results in pH values of 7.77 or 7.73, respectively. Conversely, in order to completely compensate for the CO₂-induced pH decline resulting from an atmospheric pCO₂ increase to ∼700 µatm in the A1B scenario, the surface water TA would have to increase to approximately 2800 µmol kg⁻¹ – which is a completely unrealistic TA concentration for surface water in the Gotland Sea regardless of productivity and oxygen conditions. It is for that reason beyond any doubt that the only possible way to avoid acidification of open Baltic Sea waters is to implement large reductions in CO₂ emissions. Although there is a larger pH decline in the BSAP than in the BAU scenario, the possible negative influence must be considered to be of a marginal importance compared to the vast benefits for Baltic Sea ecosystems following reduced deep water dead zones.