Declining ocean chlorophyll under unabated anthropogenic CO2 emissions

Photosynthetic assimilation of carbon dioxide and inorganic nutrients by phytoplankton constitutes a necessary prerequisite for sustaining marine life. This process is tightly linked to the concentration of chlorophyll in the ocean’s euphotic zone. According to a recent field study marine chlorophyll(a) concentrations have declined over the last century with an estimated global rate of 1.0% of the global median per year. Here we attempt to identify possible mechanisms which could explain such trends. We explore these questions using an ocean general circulation model forced with documented historic and projected future anthropogenic emissions of carbon dioxide according to the IPCC SRES A1FI emission scenario until the year 2100. We further extend the time period covered by the A1FI scenario by assuming a linear decline in emissions from 2100 to 2200 and keeping them at zero levels until 2400. Our numerical simulations reveal only weak reductions in chlorophyll(a) concentrations during the twentieth century, but project a 50% decline between 2000 and 2200. We identify a local and a remotely acting mechanism for this reduction in the North Atlantic: (I) increased sea surface temperatures reduce local deep mixing and, hence, reduce the nutrient supply from waters at intermediate depths; (II) a steady shoaling of the Atlantic overturning cell tends to transport increasingly nutrient depleted waters from the Southern Hemisphere toward the north, leading to further diminishment of nutrient supply. These results provide support for a temperature-driven decline in ocean chlorophyll(a) and productivity, but suggest that additional mechanisms need to be invoked to explain observed declines in recent decades.


Introduction
Ocean ecosystems inherently rely on the net primary production (NPP) of phytoplankton which provides the energetic base for higher trophic levels. It is estimated that nearly half of the global NPP amounting to about 45-50 Petagrams of carbon per year (Pg C yr −1 , Field et al 1998 is of marine origin, as inferred from the satellite-based Coastal Zone Color Scanner (CZCS, Feldman et al 1989) data between 1979 and 1986. Routinely, the remotely sensed ocean color data are used to derive surface chlorophyll(a) concentrations (Chl a ) from which NPP can be assessed. By comparing Chl a and NPP data inferred from the subsequently launched Sea-viewing Wide Field-of-view Sensor (SeaWiFS) satellite program spanning the years from 1997 to 2002, Gregg et al (2003) discovered a 6% decrease in global NPP since the CZCS era of the early eighties. However, by studying the evolution of the remotely sensed SeaWiFS Chl a between 1998 and 2003, Gregg et al (2005) found a global increase by more than 4%. It should be noted that this increase was mainly bound to coastal areas, while the oligotrophic subtropical gyres revealed a substantial decline in Chl a , a result also confirmed by later studies (see for example Polovina et al 2008, Irwin andOliver 2009).
In an attempt to identify the controlling mechanisms of Chl a and NPP decadal fluctuations in the ocean, Behrenfeld et al (2006) detected a strong inverse correlation between global NPP and sea surface temperature (SST) in the SeaWiFS data. A similar negative relationship was observed by Martinez et al (2009) by comparing SeaWiFS and CZCS data, and related top decadal-scale climate variability in ocean temperature and stratification. To resolve longer-term trends in Chl a NPP, and their causes, however, the length of the CZCS and SeaWiFS satellite record (1979-86 and 1997 to present) has been found insufficiently short. Boyce et al (2010Boyce et al ( , 2011 have recently provided a global compilation of annual Chl a data covering the years from 1899 until 2008. This long-term data series derived from Secchi depth measurements and in situ data is suited to the investigation of long-term climatological trends. While SSTs increased between 1950 and 2008 due to anthropogenic emissions of greenhouse gases, a drop in the globally averaged Chl a concentrations by about 1.0% of the global median per year did occur during the same time interval. Declines were particularly evident in the Atlantic and Pacific, but not the Indian Ocean. Increases in SST and related decreases in mixed layer depth and vertical mixing were discussed as primary mechanisms behind the observed declines. Another possible mechanism previously identified in modeling studies which could explain the observed Chl a decline was discussed by Schmittner (2005). This study investigates the impact of a shutdown of the Atlantic Meridional Overturning Circulation (AMOC) due to increased freshwater fluxes in the Nordic Seas during stadial/interstadial transitions. A total collapse of the AMOC would result in a more than 50% decline of the phytoplankton stocks in the North Atlantic and a 10-20% reduced global export production due to a reorganization of the global overturning circulation. However, a recent analysis of stable nitrogen isotopes preserved in corals dwelling in the Labrador Current points toward only a small decline in the strength of the AMOC during the twentieth century (Sherwood et al 2011).
The present study attempts to identify mechanisms which could address potential causes for a continuous global decline in Chl a and delivers a projection of possible future trends by employing the three-dimensional ocean-atmosphere-seaice biogeochemistry model POTSMOM-C (see section 5). More specifically, we investigate the impact of changes in the mixed layer depth and the strength of the AMOC on the Chl a concentrations. We are using the model to project a possible future scenario assuming unabated carbon dioxide emissions according to the IPCC's SRES A1FI emission path (Nakicenovic and Swart 2000).

Design of numerical experiments
In this study we employed a three-dimensional ocean general circulation model (OGCM) which is coupled to a state-of-theart marine ecosystem model (POTSMOM-C (Hofmann and Schellnhuber 2009, see section 5). Atmospheric forcing was provided from an anomaly model of the atmospheric energymoisture balance (for details see Hofmann and Maqueda 2006) which allows the calculation of deviations from the given NCEP/NCAR climatology (Kalnay et al 1996).
Anthropogenic effects arising from a climate projection under the A1FI emission pathway (Kuhlbrodt et al 2009) have been superimposed on the NCEP/NCAR climatology. The model was integrated until a steady state was reached (spin-up of 6000 yr). The steady state defines the initial condition of the calendar year 1800. Integration from 1800 to 2000 is driven by reconstructed CO 2 emission rates (Marland et al 2007). Subsequently, the model was integrated by utilizing CO 2 emission rates inferred from the IPCC SRES business-asusual emission path, A1FI, beginning in year 2000, continuing through 2100 (Nakicenovic and Swart 2000) and extending it until year 2200 by assuming a linear decline to zero, corresponding to a release of 4075 Pg C in total. The maximum increase of the global mean SST by 7.9 • C compared to preindustrial temperatures occurs around year 2150. After 2200 we kept modeled anthropogenic CO 2 emissions at zero.
It should be noted that POTSMOM-C accounts for the mineral ballast effect, which links the settling velocity of particulate organic material (POM) to its content of mineral ballast, notably calcite (CaCO 3 ) particles (Klaas and Archer 2002, Balch et al 2010, Iversen and Ploug 2010. As biogenic calcification weakens due to CO 2 -driven ocean acidification, the remineralization depth of organic material will become shallower, and vice versa (see section 3).

Evolution of chlorophyll concentrations
In steady state the atmospheric partial pressure of CO 2 ( pCO 2 ) is determined by the air-sea-gas exchange parameterized according to Wanninkhof (1992) stabilizing at 282 μatm. The simulated NPP of 43 Pg C yr −1 and export production (100 m level) of 8.6 Pg C yr −1 compare well with observations (Field et al 1998, Najjar et al 2007. The modeled circulation patterns are close to observations, such as for the AMOC with a maximum strength of 22 Sv (1 Sv = 10 6 m 3 s −1 ) in steady state. Likewise, the preindustrial distribution of diagnosed Chl a concentrations in the world ocean reproduces the wellknown pattern of high values in polar and subpolar regions and low values in the oligotrophic subtropical gyres (figure 1(a)). Chl a concentrations were derived diagnostically from phytoplankton concentrations provided by a nutrient phytoplankton zooplankton detritus (NPZD) model Schellnhuber 2009, Six andMaier-Reimer 1996) assuming a variable phytoplankton chlorophyll to carbon (Chl a :C) ratio similar to that used by Aumont et al (2003) (see section 5).
It is assumed that Chl a concentrations in the mixed layer depend on the degree of light and nutrient limitations. As mixed layer depths (MLDs) increase, phytoplankton cells will become exposed to low light levels stimulating the production of Chl a , and vice versa. The parameterization also recognizes that severely nutrient-limited cells produce less Chl a than those under nutrient replete conditions.
In our model simulations global mean SSTs have risen by about 0.55 • C between preindustrial times and 1995. More than two-thirds of this increase (0.40 • C) occurred between 1950 and 1995, which compares well with empirical observations (Reynolds et al 2002). As a consequence,  (d)) only small differences can be recognized. Notably, some regions of the North Atlantic, the Arctic and the SO reveal a decline in Chl a surface concentrations by up to 10%, while at mid-latitudes a small increase was found.
Unabated anthropogenic CO 2 emissions, as assumed in our model simulations, resulted in a global increase of average SSTs by about 6.0 • C by year 2100 leading to a strongly stratified ocean which inhibits nutrient supply to the euphotic zone through reduced vertical mixing of subsurface and intermediate depth waters. Consequently, sea surface Chl a concentrations during the summer season declined by up to 50% globally, compared to preindustrial levels (figures 1(c) and (e)). The examination of the model output regarding the Chl a -SST relationship confirms the negative correlation between these variables previously described by other authors (Behrenfeld et al 2006, Boyce et al 2010. Figure 2 shows a scatter plot of annual mean surface SST versus Chl a for different ocean basins. Between 1900 and 2100 (black circles in figure 2) there is a ubiquitous decline in Chl a concentrations with rising SSTs on a global scale. As an example (see figure 2), basin mean Chl a concentrations in the South Atlantic decline from 1.20 mg m −3 in year 1800 to 0.9 mg m −3 in 2100 while SSTs increase by about 8 • C.
After 2100, as CO 2 emissions decrease toward zero, and SSTs slightly cool, Chl a concentrations tend to recover toward baseline values in most of the ocean basins (red dots in figure 2). This recovery of phytoplankton chlorophyll toward the end of the simulation can be mostly explained by a redistribution of the nutrients in the ocean and an increase of the MLDs. As an example, the weakening and shoaling of the AMOC leads to an accumulation of phosphate (nutrient trapping) in the subsurface region of the Southern Ocean (figures 3(b) and (d)) resulting in an intensification of the nutrient supply at the sea surface. In the North Pacific the MLD starts to recover toward its initial values (figures 4(a)-(d)), resulting in a rebound of the Chl a concentrations toward preindustrial values. Also, the effect of ocean acidification on the downward export of POM contributes to this trend in the Chl a concentrations. Since the ecosystem model accounts for the effect of mineral ballast-notably calcite particles (CaCO 3 )-in vertically distributing POM (Klaas and Archer 2002), the decline in CaCO 3 export by 75% under elevated atmospheric pCO 2 levels (Hofmann and Schellnhuber 2009) leads to a shoaling of the remineralization depth in high production zones (red line in figures 3(b) and (d)) by nearly 100 m. As a result sea surface phosphate concentrations ([PO 4 ]) in most regions start to recover after year 2200.
Although the mixed layer depth increases in most regions (figures 4(b) and (c)), recovery of Chl a is not observed to the same extent in the North Atlantic and Arctic Oceans. Here, Chl a concentrations remain at low levels (figure 2), explained by a continuously decreased phosphate availability after 2100 ( figure 5(g)). Nutrient concentrations within the mixed layer of North Atlantic and Arctic Oceans are determined by a number of physical and biogeochemical processes. The supply of PO 4  According to our model results, rising SSTs owing to increasing atmospheric CO 2 levels (figures 5(a) and (b)) will negatively affect phytoplankton by two hypothesized mechanisms: first, the shoaling of the average MLD in the North Atlantic between 50 • N and 60 • N and from 160 to 110 m depth (figure 5(c)) (between 1800 and 2100) leads to a weaker nutrient supply of surface waters from intermediate depths and, hence, declining surface Chl a concentrations ( figure 5(h)). Around the middle of the twenty-first century the AMOC strength starts to decline remarkably (figure 5(d)) in concert with a shoaling of the overturning cell by about 1000 m ( figure 3(d)). The latter results in uplifting of the LNM (figure 5(e)). As a consequence, the northbound PO 4 transport at the equator declines (figures 3(c) and 5(f)), resulting in a time-retarded lowering of the North Atlantic PO 4 concentrations between 50 • N and 60 • N in the mesopelagic zone ( figure 5(g)), entailing a drop in the North Atlantic Chl a concentrations by nearly 50% on average by year 2200.

Discussion and conclusions
Our numerical simulations aimed at a deeper understanding and a quantitative assessment of the sensitivity of sea surface Chl a concentrations to rising atmospheric pCO 2 levels from year 1800 until 2400. We have restricted our discussion mostly to the analysis of the impact of ocean physics on Chl a concentrations, because the biological mechanisms are uncertain. Although the model simulations confirm the observed twentieth century trend and its relation to SSTs, the magnitude of the simulated Chl a decline is much lower than that observed by Boyce et al (2010Boyce et al ( , 2011. The reason for that observed strong decline in Chl a concentrations during the second half of the twentieth century remains unresolved by our model. Multiple possible causes, which have been ignored in our simulations, such as an increased grazing pressure by zooplankton, a change in the light limitation due to a lower cloudiness, changes in phytoplankton species composition, enhanced sea-ice melt, among others, could also affect Chl a concentrations over time. Future modeling studies may want to explore the sensitivity of long-term Chl a trends to these factors.
Although our model predicts only a small Chl a decline over the twentieth century, it does predict a major decline for the twenty-first century. As an example, Chl a concentrations in the North Pacific drop from about 0.8 mg m −3 in year 1900 to approximately 0.5 mg m −3 in 2100 (see figure 2) due to a shoaling of the MLD (see figures 4(a) and (b)), which reduces the supply of nutrients from the subsurface ocean. The subsequent return to nearly initial Chl a values is due to a reorganization of the ocean circulation and a deepening of the maximum MLD (figures 4(b) and (c)) combined with a shoaling of the remineralization depths (figures 3(b) and (c)) toward the end of the simulation. A permanent and most pronounced Chl a reduction, however, occurs in the North Atlantic and Arctic Oceans, which is of the order of nearly 50% of the global average by year 2200 compared to preindustrial levels. We identify two main mechanisms leading to a Chl a decline, although they are not sufficient to explain the large trend found by Boyce et al (2010Boyce et al ( , 2011: aside from a lowering of the MLDs until year 2100 owing to increasing SSTs combined with a freshening of the sea surface waters, the shoaling and weakening of the AMOC cell was identified as being responsible for the declining Chl a concentrations in the North Atlantic and Arctic regions in the next centuries. This implies a continuous nutrient depletion in the northbound water masses of the upper AMOC limb over the simulated time interval leading to a progressing nutrient limitation in the North Atlantic and Arctic Ocean areas. A previous study by Schmittner (2005) showed that a complete shutdown of the AMOC due to freshwater input from melting ice masses will lead to declining phytoplankton stocks in the North Atlantic by more than 50%. In his study, the AMOC's complete cessation accounts for an increase of nutrient concentrations in the subsurface and mesopelagic zone owing to the spread of nutrient rich Antarctic bottom water (AABW) toward the north displacing the nutrient depleted North Atlantic deep water (NADW), while sea surface nutrient concentrations dropped to almost zero due to high stratification (see figures 4(c) and (e) in Schmittner 2005). In contrast to the work by Schmittner (2005), we get a reduction of the nutrient concentrations instead of an increase in the mesopelagic and subsurface zone in the North Atlantic (see figures 3(b) and (d)). The reason for this difference is as follows: while in the study by Schmittner (2005) the AMOC shutdown allows for the invasion of nutrient rich AABW, in our simulations the shoaling AMOC imports more and more nutrient depleted Southern Ocean subsurface and mode waters toward the north. This is why even the more modest AMOC reduction by 40% leads to a drop in the annual mean Chl a concentrations by up to 50%-a number comparable to that by Schmittner (2005). Consequently, the marine ecosystem in the North Atlantic and Arctic Oceans is predicted to exhibit strong vulnerability against even a moderate weakening of the AMOC. A continuation of the observed decline in marine Chl a under business-as-usual emissions scenarios is therefore considered a real possibility that should be monitored empirically, and taken into account when developing future management scenarios.

POTSMOM-C
The POTSMOM-C model (Hofmann and Schellnhuber 2009) employed in the current study is an ocean-sea-ice-atmosphereenergy-moisture-balance general circulation model based on MOM-3 (Pacanowski and Griffies 1999), which was coupled to a biogeochemistry sub-model derived from the HAMOCC-3.1 code (Six and Maier-Reimer 1996). The horizontal resolution is 3.75 • × 3.75 • . In the vertical, the grid consists of 24 levels with thickness increasing from 25 m at the top to about 500 m at the bottom. POTSMOM-C benefits from a second order moments advection scheme invented by Prather (1986), which is nearly free of spurious diffusion and dispersion (Hofmann and Maqueda 2006).
The growth rates of phytoplankton have been parameterized as a function of photosynthetically available radiation (PAR), temperature, mixed layer depth and the concentrations of micro (iron) and macro (phosphate, nitrate) nutrients, where only the scarcest nutrient available is regarded to be the limiting factor (Liebig's law). The biogeochemistry model is isogeochemical, i.e. chemical weathering and exchange processes with the oceanic sediment are not included. In this study we utilize the results from anthropogenically forced (A1FI 2100, Nakicenovic and Swart 2000) climate simulation by Kuhlbrodt et al (2009), as described more in detail in Hofmann and Schellnhuber (2009). The wind stresses were kept constant, i.e. we only used forcing data from NCEP/NCAR (Kalnay et al 1996) even for the transient simulations.

Chlorophyll
The actual Chl a :C ratio in each grid cell is diagnosed by Here Chl M /C = 1/37 and Chl m /C = 1/90 are the maximum and the minimum Chl a :C ratios in phytoplankton in units of mg Chl (mg C) −1 , respectively. f (I 0 ) is the light limitation function defined by the depth integral over the euphotic zone of depth z: f (I 0 ) = I 0 αPAR(1/z) z dz exp{−kz }. I max PAR is the critical irradiance for photoadaptation (90 W m −2 ), I 0 is the solar irradiance at the sea surface in units of W m −2 , k is the light extinction coefficient (0.025 m 2 (W day) −1 ), PAR = 0.4 is the photosynthetically active radiation and a = 0.03 day −1 m 2 W −1 is the initial slope of the photosynthetic P-I curve (Aumont et al 2003). [PO 4 ] and [NO 3 ] are the surface phosphate and nitrate concentrations in μmol l −1 while K PO 4 = 0.016 μmol l −1 and K NO 3 = 0.013 μmol l −1 are the corresponding half-saturation constants for phosphate and nitrate uptake, respectively. For the sake of simplicity, we have neglected the dependence of Chl a on the dissolved iron concentrations.