Coupled physical / biogeochemical modeling including O 2-dependent processes in the Eastern Boundary Upwelling Systems : application in the Benguela

The Eastern Boundary Upwelling Systems (EBUS) contribute to one fifth of the global catches in the ocean. Often associated with Oxygen Minimum Zones (OMZs), EBUS represent key regions for the oceanic nitrogen (N) cycle. Important bioavailable N loss due to denitrification and anammox processes as well as greenhouse gas emissions (e.g, N 2O) occur also in these EBUS. However, their dynamics are currently crudely represented in global models. In the climate change context, improving our capability to properly represent these areas is crucial due to anticipated changes in the winds, productivity, and oxygen content. We developed a biogeochemical model (BioEBUS) taking into account the main processes linked with EBUS and associated OMZs. We implemented this model in a 3-D realistic coupled physical/biogeochemical configuration in the Namibian upwelling system (northern Benguela) using the high-resolution hydrodynamic ROMS model. We present here a validation using in situ and satellite data as well as diagnostic metrics and sensitivity analyses of key parameters and N2O parameterizations. The impact of parameter values on the OMZ off Namibia, on N loss, and on N 2O concentrations and emissions is detailed. The model realistically reproduces the vertical distribution and seasonal cycle of observed oxygen, nitrate, and chlorophyll a concentrations, and the rates of microbial processes (e.g, NH +4 and NO − 2 oxidation, NO−3 reduction, and anammox) as well. Based on our sensitivity analyses, biogeochemical parameter values associated with organic matter decomposition, vertical sinking, and nitrification play a key role for the low-oxygen water content, N loss, and N2O concentrations in the OMZ. Moreover, the explicit parameterization of both steps of nitrification, ammonium oxidation to nitrate with nitrite as an explicit intermediate, is necessary to improve the representation of microbial activity linked with the OMZ. The simulated minimum oxygen concentrations are driven by the poleward meridional advection of oxygen-depleted waters offshore of a 300 m isobath and by the biogeochemical activity inshore of this isobath, highlighting a spatial shift of dominant processes maintaining the minimum oxygen concentrations off Namibia. In the OMZ off Namibia, the magnitude of N 2O outgassing and of N loss is comparable. Anammox contributes to about 20 % of total N loss, an estimate lower than currently assumed (up to 50 %) for the global ocean.


Introduction
The Eastern Boundary Upwelling Systems (EBUS) (California, Humboldt, Canary, and Benguela upwelling systems) are specific areas connecting the coastal zone to the open ocean with the subtropical gyres.The eastern part of these gyres is forced by equatorward winds, especially the trade winds.These winds drive the eastern boundary currents, specific to each of these EBUS.Due to Earth's rotation, the topography, and the coastal boundary, a wind-driven offshore Ekman flow takes place in the surface layer.This flow creates an upward vertical advection of cold and nutrient-rich deep waters.The injection of nutrients in the euphotic layer plays a major role in the development of phytoplankton and zooplankton biomasses up to fish.These EBUS represent 20 % of the global catches (Fréon et al., 2009) with only 1 % of the global surface.Other important features contribute to these highly productive zones such as the wind curl (responsible of a second upward vertical advection zone offshore of the continental shelf), the coastal topography, and the poleward undercurrent with water masses enriched in nutrient and depleted in oxygen originating from the equatorial zone.Currently, these EBUS are crudely represented in the global climate atmosphere-ocean models used within the Coupled Model Inter-comparison Project 5 (CMIP5) due to their coarse resolution.A mean warm bias of 2-3 • C (difference between the mean over the different CMIP5 models and the mean of the observations) is estimated for these models (Toniazzo and Woolnough, 2013).Improving our capability to represent these areas is of crucial importance, especially in the climate change context.Indeed, these areas might experience an important change in the winds (magnitude, direction) (i.e., Garreaud and Falvey, 2008) and thus in productivity (i.e., Bakun et al., 2010).Moreover, Stramma et al. (2008Stramma et al. ( , 2010) ) also noticed a decrease of dissolved oxygen concentrations in the tropical ocean over the last 50 yr, with a higher decrease in the Atlantic ocean.Important greenhouse gas emissions (N 2 O and CO 2 ) occur in these EBUS, as shown recently in the Humboldt system (Paulmier et al., 2008).However, future trends of these emissions are unknown.Adequate tools have to be developed to address this question.We need high-resolution models (a few km) combined with a sufficiently detailed description of the biogeochemical processes, relevant for the O, N, and C biogeochemical cycles under oxic, hypoxic (O 2 < 60 mmol O 2 m −3 ), and suboxic conditions (O 2 < 25 mmol O 2 m −3 ).Regional models offer this opportunity.
In this paper we present a regional high-resolution coupled physical/biogeochemical model taking into account the pertinent O 2 -dependent processes to follow chlorophyll a, nitrate, and O 2 concentrations as well as nitrogen (N) loss and N 2 O emissions that we applied in the Benguela upwelling system.We discuss the impact of key parameter values on these key quantities.
Among the different EBUS, the Benguela upwelling system in the South Atlantic Ocean plays a special role.This EBUS presents one of the highest primary productions of all EBUS (Carr, 2002;Carr and Kearns, 2003).It is the only one bordered by two warm currents: the Angola Current, enriched in nutrients and poor in oxygen content in the northern part, (Monteiro et al., 2006;Mohrholz et al., 2008), and in the southern part the Agulhas Current with its eddies, carrying cold (warm) and enriched nutrient (poor nutrient content) waters for the cyclonic (anticyclonic) eddies (Shannon, 2006).Baroclinic instabilities and interactions with wind and bathymetry create a large variety of mesoscale and submesoscale structures (eddies, filaments, fronts, etc.) inducing a strong mixing (Penven et al., 2001;Shannon, 2006;Veitch et al., 2009).The Benguela upwelling system is divided into two subsystems mainly due to the trade wind regime.The northern part is characterized by a permanent upwelling with seasonal intensity (maximum in winter).Due to its high productivity and subsequent export production as well as the circulation, an intense remineralization occurs between 100 and 600 m depth, creating an Oxygen Minimum Zone (OMZ), especially in the Namibian upwelling system between 20 • S and 25 • S (Monteiro et al., 2006;Hutching et al., 2009).This OMZ is controlled by local productivity and stratification as well as by remote forcing associated with the oxygenlow content of the poleward Angola Current and associated poleward undercurrent (Monteiro et al., 2008(Monteiro et al., , 2011)).In this OMZ, denitrification and anammox processes occur and induce nitrogen loss (Kuypers et al., 2005;Lavik et al., 2008).During anoxic conditions, sulphur emissions can occur with their subsequent impacts (respiratory barrier for zooplankton, high fish mortality, etc.).A warming trend was observed using the Sea Surface Temperature (SST) satellite data in this EBUS (north and south) from 1997 up to 2006-2007 -a different feature as compared to the other EBUS (Rouault et al., 2007;Demarcq, 2009).
In any modeling work, in situ observations represent an essential component and cannot be dissociated from modeling efforts.In the Benguela, thanks to several initiatives from different countries (e.g., South Africa, Namibia, Germany, Norway), observations of important variables and fluxes were made as N 2 O concentrations, nitrification, denitrification, and anammox processes.Indeed, these three processes represent a loss of bioavailable (fixed) nitrogen through the production of gaseous products (N 2 O and/or N 2 ), with the potential to affect biogeochemical cycles.So a detailed description of the microbial loop has to be included in models -for example the two steps of nitrification.In our biogeochemical model, two kinds of representations are used based on current knowledge: (1) Anammox and denitrification processes generate a N loss based on Yakushev et al. (2007) work; (2) N 2 O production is estimated using the parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) which estimates the N 2 O production under oxygenated conditions and at low-oxygen levels, mimicking the N 2 O production from nitrification and denitrification processes.
In this paper we will especially address the following questions after an evaluation of the model performance.What are the key parameters of the biogeochemical model and their influence on the OMZ representation, N losses due to denitrification and anammox processes, and N 2 O concentrations and emissions to the atmosphere?What is the relative importance of the different coupled processes maintaining the OMZ in the northern part of the Benguela upwelling system?
In the first part of the paper, the biogeochemical BioE-BUS model is detailed.Then the performance of this model is shown using the "reference simulation" (with the best parameter set found in this study) for the physics and the relevant biogeochemical variables and fluxes in the context of the OMZ.Then, we present the sensitivity analysis performed on key parameter values to improve the model performance.Finally, we discuss the influence of these key parameters on important quantities such as the volume, minimum oxygen concentrations, N losses, N 2 O concentrations, and emissions before concluding.

Hydrodynamic model
In this study we used the Regional Ocean Modeling System (ROMS; Shchepetkin andMcWilliams, 2003, 2005) in its version with the 2-way nesting capability (ROMS-AGRIF; Penven et al., 2006a;Debreu et al., 2012).This hydrodynamic model is a split-explicit and free-surface model that considers the Boussinesq and hydrostatic assumptions when solving the primitive equations.The vertical discretization follows a sigma or topography-following stretched coordinate system.No explicit lateral viscosity is added in the modeled domain except in the sponge layer at open boundaries, where it increases smoothly on several grid points.Adaptative open boundary conditions combine outward radiations and nudging towards prescribed external boundary conditions (Marchesiello et al., 2001).The vertical turbulent closure is parameterized using the KPP boundary layer schemes (Large et al., 1994).
Recently, a ROMS-AGRIF nested configuration of the Southern Africa region (SAfE for Southern African Experiment) was developed by Penven et al. (2006b) and Veitch et al. (2009).It is a 2-way nested grid AGRIF configuration, consisting in a 1/4 • coarse grid covering Indian and Atlantic regions around South Africa, extending from 2.5 • W to 54.75 • E and from 4.8 • S to 46.75 • S (SAfE coarseresolution configuration (SAfE CR), hereafter referred to as SAfE "parent" domain in Fig. 1) including and "feeding" a 1/12 • finer grid covering the northern and southern Benguela upwelling system, extending from 3.9 • E to 19.8 • E and from 12.1 • S to 35.6 • S (SAfE high-resolution configuration (SAfE HR), hereafter referred to as SAfE "child" domain in Fig. 1).These two SAfE configurations (SAfE CR and HR) have been validated, in particular the SAfE HR for the Benguela upwelling system (Veitch et al., 2009).We therefore used the SAfE HR outputs to provide the initial and open boundary conditions (physical state variables) of the Namib- ian configuration developed here (see Sect. 2.3; small domain in red in Fig. 1).

Biogeochemical model: BioEBUS
The hydrodynamic ROMS model is coupled to a Biogeochemical model developed for the Eastern Boundary Upwelling Systems, thus named BioEBUS.The advectivediffusive equation determines the evolution of a biological tracer concentration C i : On the right-hand side of the equation, the first three terms represent the advection (with u the velocity vector), the horizontal diffusion (with K h the horizontal eddy diffusion coefficient), and the vertical mixing (with turbulent diffusion coefficient K z ), respectively.The last term stands for the "source-minus-sink" term (SMS) due to biological activity.BioEBUS (Fig. 2) is a nitrogen-based model derived from the N 2 P 2 Z 2 D 2 model (Koné et al., 2005) that was successfully used to simulate the first trophic levels of the Benguela ecosystem.N 2 P 2 Z 2 D 2 model takes into account the main planktonic communities and their specificities in the Benguela upwelling ecosystem.Nitrate and ammonium represent the pool of dissolved inorganic N. Phytoplankton and zooplankton are split into small (flagellates and ciliates, respectively) and large (diatoms and copepods, respectively) organisms.Detritus is also separated into small and large particulate compartments.A cumulative layer at the sediment-water interface exists in which sinking particles are stored.Within this cumulative layer the particules cannot be advected; they just accumulate on the sea floor, without further interaction with the overlying waters.In BioE-BUS, we added a Dissolved Organic Nitrogen (DON) compartment with the source terms (phytoplankton exudation, organic excretion of zooplankton, hydrolysis of detritus) and the sink terms (ammonification of DON) following Dadou et al. (2001Dadou et al. ( , 2004) ) and Huret et al. (2005) (see Eq. 8 below).Indeed, DOM is an important reservoir of OM and plays a key role in supplying nitrogen or carbon from the coastal region to the open ocean (Huret et al., 2005).The simulated DON is the semilabile one as the refractory and labile pools of DON have too long (hundreds of years) and too short (less than a day) turnover rates, respectively (Kirchman et al., 1993;Carlson and Ducklow, 1995).An equation for nitrite (see Eq. 10) was also included in order to have a more detailed description of the microbial loop: ammonification/nitrification processes under oxic conditions (see Eqs. 21,27,28), and denitrification/anammox processes under suboxic conditions (see Eqs. 22 up to 26 and 29) using the formulation of Yakushev et al. (2007).These processes are oxygen dependent, so an oxygen equation was also introduced into BioEBUS with the source term (photosynthesis), the sink terms (zooplankton respiration, bacteria remineralization) as well as the seaair O 2 fluxes following Peña et al. (2010) and Yakushev et al. (2007) (see Eq. 12 and the formulation used for sea-air fluxes Sect.2.2.6).To complete this nitrogen-based model, N 2 O was introduced using the parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012)), which allows determining the N 2 O production under oxic conditions and at low oxygen levels (see Eq. 13).
The model contains 12 state variables (Table 1).The formulation of the SMS terms for each of the biogeochemical tracers, with parameter values (Table 2), is given by the following: (2)    2006) using the parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012)).
f Vertical extrapolation of Chl a from the surface values using Morel and Berthon (1989) parameterization (roms tools;Penven et al., 2008).

Primary production
The growth rate J P i (PAR, T , N) of phytoplankton P i (i represents flagellates or diatoms) is limited by light availability for photosynthesis (PAR: photosynthetically active radiation), temperature T ( • C), and nutrients (N represents NO − 3 , NO − 2 , and NH + 4 ), aJ P i (PAR, T ) represents the phytoplankton P i growth rate limitation by PAR and temperature, using the analytical formulation of Evans and Parslow (1985): where J max P i is the maximal light-saturated growth rate, and function of temperature (Eppley, 1972): Exponential decrease of light intensity is formulated as in Koné et al. (2005): PAR 0 is the incident radiation at the surface of the ocean, k w and k chl a are the light attenuation coefficients due to water and to chlorophyll a, respectively, θ is the chlorophyll / carbon ratio, R C/N is the carbon / nitrogen Redfield ratio for phytoplankton, [P t ] represents the sum of nano-and Table 2. Parameter values of the BioEBUS biogeochemical model."Values" are the parameter values used in the reference simulation."Original values" are from the original papers, and changed for the needs of our configuration.We recall the references used by Koné et al. (2005)    Limitation by nutrients (Fasham et al., 1990) is given by a Michaelis-Menten formulation: The phytoplankton P i growth rate is limited by NO − 3 , NO − 2 , and NH + 4 .NH + 4 is preferred to NO − 3 and NO − 2 by phytoplankton.K NO3 P i and K NH4 P i are the half-saturation constants for NO − 3 + NO − 2 and NH + 4 uptakes by flagellates or diatoms.Small phytoplankton cells are more adapted to low nutrient and stratified conditions than larger ones.Halfsaturation constants for NO − 3 and NO − 2 are usually higher than those for NH + 4 , for both flagellates and diatoms (Eppley et al., 1969;Caperon andMeyer, 1972a, 1972b).Besides the more elevated surface / volume ratio, small cells have better assimilation efficiency than large cells (Nalewajko and Garside, 1983).So, the half-saturation constants are lower for flagellates than diatoms.
Note that here it is assumed that phytoplankton is not limited by phosphate (Dittmar and Birkicht, 2001;Tyrrell and Lucas, 2002) and/or silicate and/or other micronutrient like metals (e.g., Fe), which is a reasonable assumption for the Namibian upwelling system.

Grazing
The specific feeding rate of a predator Z j on a food type P i is calculated according to the following formulation (Tian et al., 2000(Tian et al., , 2001)): g max Z j is the maximum grazing rate of the predator Z j (j represents ciliates or copepods), e Z j P i is the preference of the predator Z j to the prey P i , [P i ] is the prey concentration, k Z j is the half-saturation constant of the predator Z j for ingestion, and F t is the total biomass of available food for the predator Z j (Fasham et al., 1999).

Decomposition of particulate and dissolved detritus
Decomposition of particulate or dissolved organic nitrogen (Det represents D S , D L or DON in the following equations) is formulated as follows (Yakushev et al., 2007): (1) DcDet(O 2 ) is the decomposition of Det in oxic conditions (or ammonification): The parameterization from Yakushev et al. (2007) takes into account the temperature influence and an oxygen dependence using a Michaelis-Menten formulation: where K N is the decomposition rate in oxic conditions (K NP4 for D S and D L , and K ND4 for the DON).
(2) DcDet(NO 3 ) is the decomposition of Det in suboxic conditions (denitrification): The relative consumption of NO − 3 and NO − 2 during the classic reaction (Richards, 1965) can be determined using Anderson et al. (1982): Yakushev et al. (2007) consider the suboxic decomposition of particulate (D S , D L ) or dissoled organic matter (DON) in two stages: Then, the 1st and 2nd stages of denitrification are

Anammox
As the anammox process occurs in the Benguela upwelling system, this process is taken into account using the formulation of Yakusehv et al. ( 2007):

Fluxes at the ocean-atmosphere interface
O 2 and N 2 O fluxes at the ocean-atmosphere interface are expressed using the gas transfer velocity from Wanninkhof (1992), and the Schmidt number from Keeling et al. (1998) for O 2 and from Wanninkhof (1992) for N 2 O. O 2 saturation concentrations at one atmosphere total pressure for water-saturated air are determined using Garcia and Gordon (1992).N 2 O concentrations in equilibrium with moist air at total pressure of one atmosphere come from Weiss and Price (1980).The dry mole fraction of atmospheric N 2 O is assumed to be 318 ppb (Lueker et al., 2003;Cornejo et al., 2006;Anonymous, 2008).

Chlorophyll/nitrogen ratio
Chlorophyll a concentrations ([Chl a]; in mg Chl m −3 ) are derived from simulated phytoplankton concentrations ([P]; in mmol N m −3 ) assuming a variable Chl / N ratio following (Hurtt and Armstrong, 1996): 1.59 represents the standard Chl / N ratio (in g Chl (mol N) −1 ).If growth is light limited, then the Chl / N ratio is maximum and χ = χ max = 1, hence C/Chl min = 50 g C (g Chl) −1 .If phytoplankton growth is nutrient limited, χ = nutrient-limited growth rate/lightlimited growth rate, and the upper limit for the (C/Chl) max is fixed to 160 g C (g Chl) −1 (Charria et al., 2008).In this case the applied χ ratio increases with the quantity of light available for a constant growth rate.

Sensitivity analysis
Sensitivity analyses were performed on key parameters in order to better represent the distribution of the well-known biogeochemical fields (O 2 , NO − 3 and Chl a) and the rates of the microbial loop (NH + 4 and NO − 2 oxidation, NO − 3 reduction, anammox).Different N 2 O parameterizations were also tested to study the impact on the N 2 O concentrations in the OMZ.The objective is to evaluate the sensitivity of key biogeochemical parameters on the OMZ off Namibia, on N losses due to denitrification and anammox processes and on N 2 O emissions to the atmosphere.The sensitivity analyses are detailed in Sect. 4. The most adjusted simulation represents the "reference simulation" and is confronted to data in the following section (Sect.3).

Namibian configuration
The Namibian configuration (Fig. 1) is built on a Mercator grid, spanning 5 • E to 17 • E and 19 • S to 28.5 • S with a horizontal resolution of 1/12 • (ranging from 8.15 km in the south to 8.8 km in the north).The grid has 32 sigma levels stretched so that near-surface resolution increases.In the coastal area, with a minimum depth of 75 m, the thickness of the first (to the bottom) and last (near the ocean-atmosphere interface) levels is 11.4 m and 0.4 m, respectively.In the openocean area, with a maximum depth of 5000 m, the thickness of the first and last levels is 853.5 m and 5 m, respectively.Bottom topography from 1 GEBCO (General Bathymetric Chart of the Oceans: http://www.gebco.net)(Fig. 1) product was interpolated onto the model grid and smoothed as in Penven et al. (2005) in order to reduce pressure gradient errors.The ROMSTOOLS package (Penven et al., 2008, http://www.romsagrif.org) was used to build the model grid, atmospheric forcing, and initial and boundary conditions.
The Namibian configuration used in this study is a small domain which is a subdomain of the SAfE HR climatological configuration (Veitch et al., 2009) (Fig. 1).For temperature, salinity, free surface, and the velocity (zonal and meridional components), the initial conditions come from the 1 January of the 10th year of SAfE HR (Veitch et al., 2009); the open boundary conditions are also provided by the 10th year of this simulation for which variables were averaged every 5 days.
A 1/2 • resolution QuikSCAT (Liu et al., 1998) monthly climatological wind stress (courtesy of N. Grima, LPO, Brest, France) based on data spanning from 2000 to 2007 is used to force the model at the surface.Surface heat and freshwater fluxes are provided by 1/2 • resolution COADS-derived monthly climatology (Da Silva et al., 1994).An air-sea feedback parameterization term, using the 9 km Pathfinder climatological SST (Casey and Cornillon, 1999) is added to the surface heat flux to avoid model SST drift (Barnier et al., 1995;Marchesiello et al., 2003).A similar correction scheme is used for Sea Surface Salinity (SSS) because of www.biogeosciences.net/10/3559/2013/Biogeosciences, 10, 3559-3591, 2013 the paucity of evaporation-precipitation forcing fields as in Veitch et al. (2009).
For biogeochemistry, initial and open boundary conditions for NO − 3 and O 2 concentrations are provided by the monthly climatology of the CSIRO Atlas of Regional Seas (CARS, 2006).Phytoplankton (P S and P L ) is a function of Chl a derived from SeaWiFS climatology and vertically extrapolated from the surface values using Morel and Berthon (1989) parameterization.For zooplankton (Z S and Z L ), we applied a cross-shore profile following in situ observations, with low concentrations offshore (based on Koné et al., 2005) and increasing concentrations onshore.N 2 O is a function of O 2 from the CARS database ( 2006) using the parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012)).For other biogeochemical tracers, initial and lateral boundary conditions are established using a vertical constant profile (for D S and D L ) or exponential profile (for NO − 2 , NH + 4 and DON) (see Table 1) based on Koné et al. (2005).
The simulation is run for a total of 19 yr.A physical spin-up is performed over 7 yr as the model needs a few years to reach a stable annual cycle, then the coupled physical/biogeochemical model is run for 12 yr.The last 8 yr (years 12-19) will be analyzed in the following sections of the paper.Over these 8 yr the coupled model presents balanced nitrogen and oxygen budgets.The time step of the physical model is 15 min, and the biological time step is 5 min.Outputs are saved as 3-daily averages.

Data used
In the studied area, simulated fields are compared with different types of available data sets.The CSIRO Atlas of Regional Seas allowed establishment of monthly climatology with 1/2 • of resolution for temperature and salinity fields (CARS, 2009) and oxygen and nitrate concentrations (CARS, 2006).The World Ocean Atlas (WOA) 2001 includes a seasonal and annual global climatology of Chl a (Conkright et al., 2002).The monthly climatology SeaWiFS products of level 3-binned data (9 km, version 4, O'Reilly et al., 2000), from 1997 to 2009, processed by the NASA Goddard Space Flight Center and distributed by the DAAC (Distributed Active Archive Center) (McClain et al., 1998), are used for comparison with surface simulated Chl a. Simulated Chl a is the sum of flagellate and diatom concentrations expressed in N units and converted into Chl a (mg Chl m −3 ) using a variable ratio described (see Sect. 2.2.7,Eq. 30).
Different in situ data (sections and stations) are used for the model/data comparison.Their positions are shown in Fig. 3. Temperature, salinity, oxygen, nutrients, Chl a, primary production, and mesozooplankton data were collected in May 1998 during the AMT 6 cruise (Aiken, 1998;Aiken and Bale, 2000;Aiken et al., 2000).Samples were collected during expeditions M57/2 of RV Meteor February 2003 (Zabel et al., 2003;Kuypers et al., 2005) and AHAB1 of RV Alexander von Humboldt in January 2004 (Lavik et  al., 2008), and allowed estimation of temperature, salinity, oxygen, and nutrients in the Namibian upwelling system.In October 2006 the Danish Galathea expedition crossed the Benguela upwelling system.Measurements of temperature, salinity, oxygen, nutrients, and Chl a were collected (courtesy of L. L. Søerensen, National Environmental Research Institute, Denmark) at different stations and along a vertical section.The AMT 17 cruise in November 2005 estimated DON concentrations in the southern Benguela.For the first time in the Namibian upwelling system, N 2 O samples were collected off Walvis Bay at around 23 • S during the FRS Africana cruise in December 2009 within the framework of the GENUS project (Geochemistry and Ecology of the Namibian upwelling system).The collected samples were poisoned with mercuric chloride onboard and measurements were done after the cruise at IFM-GEOMAR, Germany, by using the static equilibration method according to Walter et al. (2006).
All data enumerated above will be compared to simulations outputs for the same geographical positions (except for DON because there are no available measurements in the Namibian upwelling system) and the same climatological months.

Model evaluation
To compare simulated fields m of the "reference simulation" and data d (in situ, satellite, and climatological data), different statistical metrics are selected: the mean (M), the bias (M m − M d ), the root mean square (RMS), the standard deviation (σ ), and the correlation coefficient between simulated fields and data (R) that are combined in a Taylor diagram (Taylor, 2001).For the statistics the model is interpolated onto the observed data locations.The performance of the coupled model is studied over the whole 3-D simulated domain (Namibian configuration).The model is in agreement with data if R ≥ 0.9, 0.6 ≤ σ * ≤ 1.4 and RMS * ≤ 0.4 (R: correlation coefficient, σ * : normalized standard deviation, and RMS * : normalized centered pattern RMS difference).

Temperature, salinity, and density
We first evaluate the model capability to represent temperature, salinity and density fields in the Namibian upwelling system.Density is computed from temperature and salinity data using the same Jackett and McDougall (1995) relationship as the one used in ROMS.The statistics are computed between the simulated fields and in situ data as well as annual and seasonal CARS (2009) climatologies.For comparison with seasonal CARS (2009) climatology, statistics are performed between the surface and 600 m depth as below this depth they do not show clear seasonal variation.The annual fields are also compared in the top 600 m as well as over the whole water column.
For the three fields, statistics highlight an important variability when model outputs are compared with in situ data.The correlation coefficient usually stays above 0.9, the normalized standard deviation is between 0.6 and 1.4, and the normalized centered pattern RMS difference is less than 0.45 (Fig. 4a, b, and c).In general, the model outputs slightly overestimate the temperature field, and underestimate salinity and density fields as compared to in situ data for specific years (Table 3).The most important differences between simulated and in situ data are mainly observed in the surface layer (0-100 m) over the continental shelf and slope because this area is very sensitive to the different forcing fields.However, no clear latitudinal trend can be highlighted.These results especially point out an important interannual variability for the in situ data that is not represented using our simulation.Indeed, the monthly climatological fields (winds, heat, and freshwater fluxes) used to force the model are linearly interpolated between each month and are repeated each year, so the forcing fields are smoothed in time.
As the modeled outputs represent a climatological situation, and not the interannual variability, the same statistics are now presented between the simulated fields and CARS ( 2009) climatologies (green and red symbols in Fig. 4).The correlation coefficient is above 0.96 for three fields.The normalized standard deviation is between 0.93 and 1 for temperature, between 0.83 and 0.93 for salinity, and between 0.9 and 1 for density.The normalized centered pattern RMS difference is less than 0.2 (or a RMS difference of 0.78 • C) for temperature, 0.27 (or a RMS difference of 0.096) for salinity, and 0.28 (or a RMS difference of 0.16 kg m −3 ) for density.The biases present the same trends as those observed with in situ data: a slight overestimation of temperature in our simulation (positive biases between 0.053 • C during spring and 0.33 • C during winter), while salinity and density fields are slightly underestimated (negative biases between −0.043 during spring and −0.0063 during winter, −0.069 kg m −3 during winter and −0.017 kg m −3 in autumn, respectively) (Table 3).These different statistical metrics giving satisfying results, we illustrate in the following paragraph two types of model/data distributions: in space (Fig. 5) and in time (Fig. 6).
The simulated 8 yr mean of temperature, salinity, and density is compared to the annual climatology from CARS database (2009) along a vertical section (west-east direction) at 23 • S (Fig. 5).As a signature of the coastal upwelling off Namibia, simulated temperature, salinity, and density isocontours rise to the surface over the continental shelf and slope, in agreement with the CARS climatology (2009).In the top 100 m near the coast, simulated salinity is lower than the climatological data.This bias seems to come from the SSS corrections used in the model configuration.Indeed, the SSS values from the COADS climatology used for the correction scheme in the freshwater fluxes (Sect.2.3) are lower  .The time series (Fig. 6) confront the seasonal cycle of simulated temperature and salinity with CARS monthly climatology (2009).Both data sets follow the same seasonal cycle, with an isotherm and isohaline rise during the winter season associated with the seasonality of the trade winds in this area.The SST varies between 12-13 • C in austral winter and 17-18 • C in austral summer.The isotherm 9 • C is located around 300-350 m depth in the simulated and climatological fields.For salinity, a subsurface maximum salinity observed during summer (in the first 200 m depth) rises to the surface in austral winter, progressively, in agreement with the climatology.This subsurface maximum comes from strong poleward undercurrent along the Namibian Shelf in austral summer that feeds the Namibian system in salty Angola water masses (Monteiro and van der Plas, 2006).When the trade winds become stronger during austral winter, the dominant flow on the Namibian Shelf is equatorward (Monteiro and van der Plas, 2006) and these salty subsurface waters upwell to the surface.This feature is also present in the 10 yr time series off Walvis Bay (see Fig. 3 for location) at 23 • S (Monteiro and van der Plas, 2006).An underestimation of salinity over the first 100 m depth can be noticed for the simulated field as explained previously.Despite a lower standard deviation as compared to CARS database (2009), these statistics show that the model is able to simulate the annual mean and seasonal cycle for temperature, salinity, and density fields over the Namibian upwelling system.

Oxygen and nutrients
Simulated oxygen and nitrate fields are first compared to in situ data.As for the previous evaluated fields (Sect.3.1), the different statistical metrics highlight an important variability between the in situ data and the modeled fields.For oxygen and nitrate, the correlation coefficient is above 0.8 for all comparisons, except for AHAB1 data in January 2004 (Fig. 4d and e).The normalized centered pattern RMS difference is usually lower than 0.6 and the normalized standard deviation stands between 0.4 and 1.1, except for AHAB1 data.This standard deviation stays lower than 1 for the oxygen concentrations.The bias between the simulated and in situ oxygen concentrations is usually positive, with higher bias for AHAB1 cruise up to 67.7 mmol O 2 m −3 (Table 3).This positive bias is mainly due to the overestimation of oxygen concentrations on the continental shelf especially for water depth below 130 m.This explains the lower variation of the simulated oxygen concentrations as compared to in situ measurements (σ * < 1).From these results it appears that the model does not succeed in representing extremely low nitrate and oxygen concentrations on the bottom shelf waters.Indeed, sampling conducted during the AHAB1 cruise reveals particular conditions (important sulfidic event) in the continental shelf waters off Namibia (Lavik et al., 2008).Sulfate reduction produces toxic hydrogen sulfide under oxygen-and nitrate-depleted conditions, usually linked with the continental shelf sediment processes, which are not considered in the BioEBUS model at the moment.Except these extreme events, the seasonal cycle and annual mean of simulated oxygen and nitrate concentrations give satisfying results as compared to CARS database (2006) (Fig. 4d and e).The correlation coefficient is above 0.9 and 0.95 for oxygen and nitrates, respectively, for all comparisons with CARS.The normalized standard de-viation is between 0.84 and 0.89 for oxygen, and between 0.92 and 1.04 for nitrates.The biases vary between −5.56 and 2.95 mmol O 2 m −3 for oxygen, and between −1.89 mmol N m −3 and 1.01 mmol N m −3 for nitrates and are found in austral winter and summer, respectively (Table 3).The vertical section at 23 • S (Fig. 5d and e) and the time series (Fig. 6) confirm the good statistics for oxygen and nitrate fields as compared to CARS database (2006).The model simulates a realistic vertical structure for oxygen concentrations, with the OMZ located along the continental slope (Fig. 5d).More oxygen-depleted waters are present in the northern part of the domain (between 19 • S and 24-25 • S) where the South Atlantic Central Waters (SACW; Mohrholz et al., 2008) are present.For nitrates the simulated and climatological surface concentrations are lower than 1 mmol N m −3 in the open ocean and higher than 16 mmol N m −3 along the coast.The pool of high nitrate concentrations (around 30-35 mmol N m −3 ) is located at depths between 200 and 1500 m (with variations depending on latitude and related to the upwelling cells along the coast; not shown).Discrepancy is visible with too high oxygen and nitrate concentrations simulated over the continental shelf, and especially for water depth below 130 m, probably due to biogeochemical sediment processes, as explained before.
The time series (Fig. 6) confront the seasonal cycle of simulated oxygen and nitrate concentrations with CARS monthly climatology (2006).Both data sets follow the same seasonal cycle, with shallow deoxygenated waters in austral summer and a thicker oxygenated surface layer during austral winter.This feature is due to the alternative dominance of the poleward undercurrent in austral summer and the coastal Benguela Current in austral winter, as explained in the previous section.However, the OMZ is more expanded (vertically and temporally) in the CARS climatology (2006) than in the simulated field, with a thicker oxygenated surface layer in the latter.This difference comes from the very coastal deoxygenated water column (see above), more developed in the observational data (CARS, 2006) than in the model, in which the seasonal variations are lower Table 4. Comparison between simulated and in situ data in the modeled domain for nitrite and ammonium concentrations (mmol N m −3 ), total primary production (g C m −2 d −1 ) integrated over the euphotic zone, annual primary production (10 2 g C m −2 yr −1 ), and mesozooplankton (mmol N m −2 ) integrated over 200 m depth., 2006).The simulated oxygen concentrations are between 225 and 270 mmol O 2 m −3 close to the surface, and frequently lower than 50 mmol O 2 m −3 at depth, in agreement with this in situ time series at 23 • S.However, concentrations below 25 mmol O 2 m −3 are scarcely obtained at this location using the coupled model, while they are frequently observed in the in situ time series.An upwelling of nutrient-rich waters occurs during the austral winter season (Fig. 6d and h) linked with the wind-driven upwelling.However, a seasonal cycle is not so marked in the simulated nitrate due to important nonlinear processes associated with eddy activity simulated using the eddy-resolving ROMS model (Gutknecht et al., 2013) as compared to coarse resolution climatology.
The other nutrients, i.e., nitrite and ammonium, present low values in the simulation, with a maximum value of 0.37 and 0.56 mmol N m −3 , respectively, with a subsurface maximum usually observed around 50 m depth above the continental shelf (not shown).These two nutrients are quickly transformed into nitrate by nitrification processes.As reported in Table 4, simulated values are within the range of in situ concentrations reported in the Namibian upwelling system.During the AHAB1 cruise, both nutrients reached very high concentrations as compared to other cruises, up to 3.4 and 8 mmol N m −3 for nitrite and ammonium concentrations, respectively.Note that nitrite and ammonium concentrations respectively above 1 and 4 mmol N m −3 only correspond to extremely low oxygen concentrations (O 2 < 4.5 mmol O 2 m −3 ) on the bottom shelf waters (water depth usually below 130 m).As explained above, these extreme conditions are not represented by our climatological simulation.More in situ measurements of nitrite and ammonium concentrations are needed to conclude.However, simulated fields seem to have the same order of magnitude as in situ data, possibly at the lower limit of the measured concentrations.
The comparison between simulated fields and the CARS database ( 2006) confirms that our model correctly represents the annual mean vertical distribution and seasonal cycle of oxygen and nitrate concentrations.Vertical gradients of both fields as well as horizontal distributions give satisfying results, except concentrations near the sediment-water interface over the continental shelf for water depth lower than 130 m.The simulated nitrite and ammonium concentrations are within the range of available in situ data.

Chlorophyll a and primary production
The model is able to reproduce the low and high chlorophyll a concentrations (Chl a) offshore and along the coast, respectively, as well as a tighter onshore-offshore gradient in the southern part of the domain, and slacker one in the northern part (Fig. 7a and b).However, the amplitude between the extremes of Chl a (offshore and coastal values) is lower in the modeled field, explaining the normalized standard deviation lower than 1 in Fig. 7e.Indeed, simulated surface Chl a are around 0.3 mg Chl m −3 offshore, being slightly higher than those derived from SeaWiFS sensor (0.2 mg Chl m −3 ), and Chl a increases up to 7 mg Chl m −3 at the coast in the model, while coastal concentrations reach 10 mg Chl m −3 in the satellite data (Fig. 7a and b).The time series (Fig. 7c and 7d) highlight the good agreement between the seasonal cycles of simulated Chl a and of the SeaWiFS climatology.Between 22-24 • S, minimum Chl a are simulated in austral winter during the maximum upwelling season.Some weeks after the maximum upwelling, Chl a increases up to maximum values in austral summer.
The Taylor diagram presented in Fig. 7e gives a quantitative comparison between the simulated and observed (climatological and satellite data) Chl a patterns.The statistics are better when the simulated Chl a is compared with the more recent SeaWiFS data (1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009) than surface WOA (2001) data, being due to the biomass temporal changes observed in EBUS (Demarcq, 2009).Indeed, the northern part of the Benguela upwelling system presents a positive anomaly of surface Chl a about +0.19 mg Chl m −3 for the last 10 yr (Demarcq, 2009).The annual and seasonal mean fields present a correlation coefficient between 0.8 and 0.9 compared with SeaWiFS climatology, corresponding to the same decade as the wind used to force the model.The normalized standard deviation lower than 1 highlights an underestimation of the simulated Chl a variance in an annual, seasonal, and spatial way, and the normalized centered pattern RMS difference is less than 0.6 (or a RMS difference less than 0.22 mg Chl m −3 ).Satellite Chl a seasonal cycle is well depicted by the model (Fig. 7c, d), with statistical metrics for seasonal distribution close to those obtained for the annual mean pattern (Fig. 7e).
To complete our model/data comparison, the simulated distribution of Chl a is now compared with in situ data from the AMT 6 cruise in May 1998, where vertical profiles of Chl a were estimated from four different sources (E.Fernandez, P. Holligan, R. Barlow, and the British Oceanographic Data Center; Barlow et al., 2002Barlow et al., , 2004) at three stations situated in the middle of our domain (St.14, 15, and 16; Fig. 8).Vertical profiles of simulated Chl-a are within the range of the different estimations, with a subsurface maximum about 1.7 mg Chl m −3 observed in the first 50 m depth of the wa- ter column.The vertical pattern as well as the range of the simulated Chl a is quite in agreement with in situ data.
The simulated total primary production is now compared with the BENEFIT data from Barlow et al. (2009) (Table 4).In our modeling experiment, following the seasonal cycle of Chl a (see above), the Namibian system appears to be twice as productive in austral summer (up to 2.3 g C m −2 d −1 near the coast) than in austral winter (up to 1.2 g C m −2 d −1 ).Along the Namibian coast, Barlow et al. (2009) describe this same seasonal cycle, but the total primary production is much higher: between 0.39 and 8.83 g C m −2 d −1 during the period of February-March 2002, and between 0.14 and 2.26 g C m −2 d −1 during the period of June-July 1999, with important spatial variations.However, the 8 yr mean of simulated primary production is in agreement with satellite estimations (Table 4).In the Benguela upwelling system, annual primary production was estimated between 3.7 × 10 2 and 9.5 × 10 2 g C m −2 yr −1 (Ware, 1992;Carr, 2002;Tilstone et al., 2009, see Table 4).Brown et al. (1991) estimated the primary production in the northern Benguela system at 4.3 × 10 2 g C m −2 yr −1 inshore of the 500 m isobath.Taking this same limit (500 m isobath), the simulated estimation is equal to 5 × 10 2 g C m −2 yr −1 , in the range of previous estimations.
The spatio-temporal variability of Chl a and primary production is quite well captured in our modeling experiment.

Zooplankton and DON
The temporal evolution of simulated copepod (mesozooplankton: large zooplankton or Z L in the model) biomass at 23 • S is compared with in situ data from  Kreiner and Ayon (2008) (Table 4).To convert the copepod abundance (number of individuals per m 2 : no m −2 ) from Kreiner and Ayon (2008) in biomass (mmol N m −2 ), we used a mean dry weight of individual animal of 40.75 µg for each individual.It corresponds to a mean dry mass for Calanoides carinatus (juvenile and adult stages) as estimated by Huggett et al. (2009), a major copepod species in the Namibian system.The abundance was converted into terms of carbon by assuming that 40 % of the dry weight was carbon (Peterson et al., 1990).The C / N ratio of 5.36 used for mesozooplankton biomass comes from AMT 6 estimations in this area (Aiken, 1998;Aiken and Bale, 2000;Aiken et al., 2000).Integrated over 200 m depth, in situ data between 2000 and 2007 present an important interannual variability, with biomass from 21.5 to 150 mmol N m −2 .However, the maximum biomass is usually observed during the first part of the year (summer-early autumn) and is located between 10 (14.2 • E) and 50 nautical miles (13.5 • E) off the coast (Kreiner and Ayon, 2008).For the same transect (23 • S), our coupled model simulates a copepod biomass varying from 21.5 to 120 mmol N m −2 , with a maximum biomass during summer-early autumn, situated at the same distance from the coast as in situ data.Seasonal cycle of zooplankton follows that of Chl a with a lag time of 3 to 8 weeks necessary for zooplankton to react to the phytoplanktonic biomass increase (Postel et al., 1995).
During AMT 6 cruise, zooplankton biomass was estimated at three stations (St.15, 12, and 7; see Fig. 3 for location).Data from different size classes (200-500 µm, 500-1000 µm, 1000-2000 µm, > 2000 µm) were integrated over the first 200 m depth of the water column.To compare with our simulation (Table 4), we only take into account the fraction between 200 and 2000 µm, usually considered as mesozooplankton (copepods in the model).For the same climatological month (May), the monthly mean concentrations are within the range of in situ data but do not reach the high in situ values; they still differ by a factor of 2. However, using the 3-day average outputs for May, the simulated copepod concentrations can reach 107 mmol N m −2 , closer to the maximum concentration measured during AMT 6 cruise.Regarding the full analyzed period and domain, mesozooplankton regularly reach concentrations of about 130-145 mmol N m −2 at the beginning of the year.So, simulated copepod distribution seems to corroborate in situ measurements from Kreiner and Ayon (2008) and AMT 6 cruise, with slightly lower concentrations.
In situ measurements of DON in the Namibian upwelling system are not available, so the simulated semi-labile DON is compared with data from the AMT 17 cruise that crossed the offshore region of the southern Benguela upwelling system in November 2005 (Holligan, 2005).To estimate the semilabile part of in situ DON measurements, we assumed the deepest value for each station as refractory DON, and we removed this refractory part (assumed constant with depth) of the vertical profile of DON.For the same climatological month, simulated vertical profiles of semi-labile DON concentrations were taken in the offshore region of the Namibian upwelling system.The simulated and in situ profiles are not directly comparable; however, in situ concentrations give an indication of the vertical distribution we can expect offshore of our domain.The simulated semi-labile DON concentrations (0.5-1.5 mmol N m −3 at the surface) are within the range of in situ concentrations (0.7-3.4 mmol N m −3 at the surface).

N 2 O distribution
During the FRS Africana cruise in December 2009, intense alongshore winds were observed off Namibia, leading to upwelling of nutrient-rich and oxygen-depleted subsurface waters.This event was intense enough to make the OMZ intersect the surface along the coast around Walvis Bay (Mohrholz et al., 2009 in the cruise report from Verheye and Ekau, 2009).For the climatological month of December, our coupled model also simulates an upwelling along the coast of Namibia.
The in situ N 2 O concentrations reach up to approximately 40 × 10 −3 mmol N 2 O m −3 in oxygen-depleted wa-ters in the water column onto the shelf, and at the shelf break near the water-sediment interface and tend to decrease with increasing O 2 concentrations (Fig. 9a and b).The simulated N 2 O concentrations, averaged for climatological December month, are compared to these first data available off Namibia.The simulated fields and in situ data do not agree well because we compare climatological fields with in situ measurements for a specific year (2009).However, here we want to show that the simulated N 2 O concentrations have similar values as compared to the first in situ N 2 O data when the simulated O 2 concentrations are close to the measured O 2 concentrations.This is the case for the oxygenated water column and for waters close to the sediment onto the continental slope, with simulated values up to 30 × 10 −3 mmol N 2 O m −3 and in situ values up to 40 × 10 −3 mmol N 2 O m −3 .The simulated 3-day averaged N 2 O concentrations follow the same trend as in situ measurements for the full analyzed period (Fig. 9c).Simulated values reach 90-100 × 10 −3 mmol N 2 O m −3 .However, more in situ measurements are needed to validate our fields at low oxygen concentrations.

Sensitivity analyses of key parameters
Due to the different modifications and additions made from the original biogeochemical model of Koné et al. (2005) (see Sect. 2.2), we had to adapt the BioEBUS model, especially some of the parameter values within the range of values found in the literature (Table 2).A first adjustment was made using BioEBUS in zero spatial dimensions (time dependent only) as the coupled high-resolution model requires a lot of computer time for a simulation.For the adjustments we used the main features for the different state variables in the coastal area and in the open ocean, deduced from the data (Sect.2.4) or from the literature.As a more detailed representation of the nitrification was introduced in BioEBUS as compared to Koné et al. (2005), some of the parameter values associated to phytoplankton nutrient uptake (half-saturation constants K NO3 P S and K NH4 P L ) were adjusted in order to obtain better agreement for the Chl a and nitrate concentrations between the simulated field and data.The maximum growth rate of the large phytoplankton (diatoms) had to be reduced from 0.8356 d −1 (Koné et al., 2005) to 0.6 d −1 (Fasham et al., 1990;Oschlies and Garc ¸on, 1999;Huret et al., 2005) to obtain a better agreement between the simulated primary production and the BENEFIT data from Barlow et al. (2009) as well as the data from the AMT 6 cruise, especially in the coastal area where diatoms represent the dominant species.We also changed the grazing function with the introduction of the preference parameter for the zooplankton for different types of plankton using the formulation of Dadou et al. (2001Dadou et al. ( , 2004)), the associated maximum growth rates and half-saturation constants for ingestion, and the small zooplankton specific excretion rate.Finally, parameters linked to the microbial loop (hydrolysis of detritus µ D S and µ D L ) were tuned in order to represent realistic mean value of DOM.This first adjusted version of the model was tested in the 3-D coupled configuration, and these results represent the "initial simulation" (or Test no. 1 in Table 5; the solid black line in Fig. 10).
After this first coarse adjustment, sensitivity analyses were performed on key parameters using the 3-D coupled model in order to improve the distribution of oxygen, nitrate, and Chl a concentrations, as well as key processes of the microbial loop when estimations based on in situ measurements are available.Moreover, different N 2 O parameteriza-  6).
tions were also tested to improve the N 2 O concentrations in the deep ocean.The objective is also to evaluate the sensitivity of biogeochemical parameter values on the OMZ off Namibia, N losses due to denitrification and anammox processes, and N 2 O emissions to the atmosphere.The methodology consists in arbitrarily changing the value of key parameters from the "initial simulation" (Test no. 1).Parameter values are changed independently and their reference values are increased or decreased by a factor of 2, 4, 6, 10, and 30, and then are tested in combined sets using the same approach as Charria et al. (2008).Relevant simulations are summarized in Table 5.The coefficient values before and after adjustments are mentioned in Table 2.For each analysis, the perturbed simulated fields were compared with the annual climatologies of nitrate and oxygen concentrations from CARS (2006)    Final simulations (combination of previous simulations; see parameters values above) 14 (= Tests no. 13 4 and NO − 2 oxidation rates, NO − 3 reduction, and anammox rates were compared to in situ estimations (Fig. 11).We keep the parameter value when the difference between the simulated fields and the data is the smallest among the different tests we made, with realistic parameter values.These sensitivity analyses are illustrated in Fig. 10 and 11, in which simulated fields or rates are the 8 yr mean, spatially averaged between 22-24 • S and 10-15 • E. For N 2 O concentrations, different parameterizations were tested (Nevison et al., 2003;Suntharalingam et al., 2000Suntharalingam et al., , 2012;; see Fig. 12) before being included in the two last tested combined simulations (Tests no.14 and 15).

Principal simulated variables
The adjustments made on decomposition rate of particulate organic matter (K NP4 ) and dissolved organic matter (K ND4 ) from Yakushev et al. (2007) formulations and sinking velocity of large detritus (w D L ) clearly improve the distribution of nitrate, oxygen, and Chl a concentrations, and more especially the associated vertical gradient (Fig. 10).The oxygen concentrations are higher than climatological data as the model does not simulate extremely low oxygen concentrations onto the continental shelf (see Sect. 3.2).Excluding this area (vertical profiles overlap on Fig. 10a), the simulated concentrations get closer to CARS data.By intensifying the vertical sedimentation of large detritus (Test no. 2), the maximum of decomposition rate occurs deeper in the water column.It decreases the nitrate concentrations in the surface layer (0 to 50 m depth) with improvement of the Chl a profile and the oxygen concentrations below 50 m depth.The simulated oxygen, nitrate, and Chl a vertical profiles are in better agreement with the data when combining the increased vertical velocity with the decreased decomposition rates of detritus and DON (Test no.4).As this combined simulation (standard simulation) gives satisfying results, we are now interested in improving the relevant rates for the nitrogen cycle (especially nitrification, denitrification, and anammox processes) without changing the distribution of oxygen, nitrate, and Chl a concentrations.In this objective a series of simulations is performed and is detailed in the following section.The relevant improvements are combined in the "reference simulation" (Test no.15).

NH +
4 and NO − 2 oxidation rates, NO − 3 reduction, and anammox rates Using the "initial simulation", fast and shallow NH + 4 regeneration leads to low nitrification rates in the OMZ (order of magnitude of 10 −4 mmol N m −3 d −1 ) as NH + 4 and NO − 2 are oxidized in the surface layer.When the oxic decomposition processes occur deeper in the water column (Test no.4: standard simulation), nitrifying activity shifts vertically, and NH + 4 and NO − 2 oxidation rates are increased by a factor of 65 and 112, respectively, in the OMZ between the "initial" and "standard simulations", while the rates for the whole water  (Nevison et al., 2003;Suntharalingam et al., 2000Suntharalingam et al., , 2012)), β coefficient (Suntharalingam et al., 2000(Suntharalingam et al., , 2012)), and sensitivity analyses on key parameters.Vertical profiles were 8 yr averaged.Vertical profiles of averaged N 2 O concentrations between 22-24 • S and 12.5-15 • E, zonal extension of the OMZ at this latitude range, overlap in (a).Red stars represent in situ measurements from the FRS Africana cruise in December 2009.Table 6.In situ and simulated estimations of NH + 4 and NO − 2 oxidation, NO − 3 reduction, and anammox (10 −2 mmol N m −3 d −1 ) in the OMZ off Namibia.Simulated estimations are 8 yr mean values from the "reference simulation" (Test no.15), computed between 22-24 • S and 10-15 column are only decreased by a factor of 1.4 for both stages of nitrification (Fig. 11).The combined set of parameters tested in the fourth simulation also activates suboxic processes nonexistent in the "initial simulation" (Fig. 11) due to improvements in the oxygen vertical distribution between Tests no. 1 and 4. For the "standard simulation", we tested nitrification, denitrification and anammox rate-limited reactions from Yakushev et al. (2007) that used a 1-D hydrophysical/biogeochemical model adapted for suboxic water column conditions (Black and Baltic seas).The anammox rate tested independently has a poor impact on microbial fluxes.However, its increase (by a factor of 10) becomes relevant when combining it with a reduction of the nitrification rate by the same factor (Test no.14; Fig. 11) as here NO − 2 is the lim-iting factor for anammox bacteria.So, a slower nitrifying activity improves anammox fluxes off Namibia, with simulated anammox flux of 1.4 × 10 −2 mmol N m −3 d −1 (Test no.14; Fig. 11) in the lower limit of in situ estimations (Table 6), but it also affects the N 2 O production linked with the OMZ (see following section).Thus, we looked for an indirect way to increase the NO − 2 pool in the OMZ and thus anammox fluxes without changing the nitrification rates.An increase of the NO − 3 reduction rate represents the other way to increase the NO − 2 pool in the OMZ (Test no.9).Combining the relevant tested parameters, the "reference simulation" (Test no.15; Fig. 11) gives satisfying results for nitrification, with 8 yr mean fluxes of 3 and 7.7 × 10 −2 mmol N m −3 d −1 for NH + 4 and NO − 2 oxidations within the simulated OMZ off Namibia, respectively, and maximum fluxes of 17.6 and 37.7 × 10 −2 mmol N m −3 d −1 , respectively, as compared to in situ estimations (Table 6) reported by Kalvelage et al. (2011) andFüssel et al. (2011).NO − 2 oxidation rate exceeds NH + 4 oxidation rate due to important input by anaerobic NO − 3 reduction to NO − 2 (mean and maximum values of 5.6 × 10 −2 and 37.3 × 10 −2 mmol N m −3 d −1 , respectively, in agreement with in situ measurements; Table 6), also discussed in Lam et al. (2009) andFüssel et al. (2011).The 8 yr mean anammox flux is low (0.1 × 10 −2 mmol N m −3 d −1 ), but note that maximum values can reach 4.3 × 10 −2 mmol N m −3 d −1 , in the range of in situ estimations off Namibia (Table 6).

N 2 O parameterization
Similarly to the sensitivity analyses performed on key biogeochemical parameters, different N 2 O parameterizations were implemented in BioEBUS (Table 5) from the "standard simulation".Nevison et al. (2003) determined the N 2 O production as a function of depth and oxygen amount consumed by nitrification.If this parameterization represents a reasonable assumption for the global ocean, it is not satisfying for EBUS and associated OMZs as they do not consider the N 2 O production due to denitrification process.The parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) relates the N 2 O production to organic matter decomposition using the oxygen consumption rate in the oxygenated ocean and at lowoxygen levels.At global scale the simulation from Suntharalingam et al. ( 2012) using β = 0.01 mol N 2 O mol N −1 and the formulation from Nevison et al. (2003) give satisfying results in terms of N 2 O concentrations and vertical distribution as compared to the compiled database (Suntharalingam et al., 2012).Nitrification is the main process affecting the marine N 2 O cycle at global scale, while denitrification can play an important role at local scale.For this reason the β coefficient accounting for the N 2 O source in suboxic waters will also be tested here as in Suntharalingam et al. (2012).Using β = 0.02 to 0.04 mol N 2 O mol N −1 , results represent an upper limit for the global N 2 O source -still reasonable as compared to the recent estimates (Suntharalingam et al., 2012).These results can then be more realistic under specific biogeochemical conditions such as those found in EBUS and associated OMZs.
Between the parameterization of Nevison et al. (2003) and the one from Suntharalingam et al. (2012) using β = 0.03 mol N 2 O mol N −1 (Fig. 12a), the mean increase of N 2 O concentrations takes place in the intermediate layer (100 to 400 m depth) for the whole analyzed domain (22-24 • S and 10-15 • E).However, improvements are highlighted when zooming in on the continental shelf and shelf break where the OMZ is developed (vertical profile overlap on Fig. 12a).The parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) using β = 0.03 mol N 2 O mol N −1 gives a significant N 2 O source when using the parameter set of the "reference simulation", providing N 2 O concentrations to the range of the few in situ observations available in the Namibian upwelling system (see Sect. 3.5 and Fig. 12b).

Discussion: implication for the OMZ and associated nitrogen fluxes
Improving the distribution of oxygen, nitrate, and Chl a concentrations, as well as key processes of the microbial loop (see Sect. 4), the best fitted simulation (reference simulation) gives satisfying results compared to the climatologies and the in situ measurements (see Sect. 3).Moreover, the N 2 O parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012)), with higher β coefficient, also improves the N 2 O concentrations in the interior ocean.We now discuss the impact of parameter values on the OMZ off Namibia, N losses due to denitrification and anammox processes, and N 2 O concentrations and emissions to the atmosphere.

OMZ off Namibia
The oceanic volume of oxygen-depleted waters is very sensitive to biogeochemical parameters used in the coupled model associated with organic matter decomposition (K ND4 and K NP4 ), vertical sinking of large particles (w D L ) and nitrification (K N42 and K N23 ) (Fig. 13).The standard deviation of the hypoxic and suboxic volumes over the 8 yr of simulation is sizeable due to the intrinsic variability of the ocean using an eddy-resolving coupled model.We first examine the hypoxic water volume using the common O 2 threshold of 60 mmol O 2 m −3 (Hofmann et al., 2011).This hypoxic water volume is very sensitive to the parameter values associated with organic matter decomposition and vertical sinking of large particles (see Tests no. 1 to 4, and no. 10 to 13 in Fig. 13).For example, its volume is multiplied by a factor of 1.8 between the "initial" and "standard simulations".Improving the distribution of oxygen, nitrate, and Chl a concentration, as well as key processes of the microbial loop (see Sect. 4), we finally kept w D L × 8, K ND4 /15 and K NP4 /3 (between Tests no. 1 and 13); these improvements allowed multiplication of the hypoxic water volume by a factor of 2.2 between the beginning and the end of the sensitivity analyses.Biogeochemical parameters for nitrification, anammox, and denitrification processes, when tested independently (see Tests no. 5 to 9) or in combination with previous improvements (see Tests no. 13 to 15), do not impact the hypoxic water volume.Finally, the best-fitted simulation (reference simulation) generates a hypoxic water volume of 2.21 × 10 3 km 3 very close to CARS database (2.32 × 10 3 km 3 ) for the same O 2 threshold.
The tested biogeochemical parameter values have the same impact on the suboxic water volume or OMZ volume (O 2 < 25 mmol O 2 m −3 ) as the hypoxic water volume -however with a higher intensity.For example, the OMZ was almost inexistent when starting the sensitivity analyses (Test 1 in Fig. 13) and became noticeable in the "standard simulation" (0.8 × 10 2 km 3 ; see Test no. 4 in Fig. 13).The OMZ volume is still doubled at the end of sensitivity analyses (1.7 × 10 2 km 3 ; see Test no. 15 in Fig. 13) but differs by a factor of 2.3 from the CARS database (3.9 × 10 2 km 3 ).However, this volume difference comes from the area shoreward of 130 m isobath.Excluding this area both estimations are very close (1.7 and 1.8 × 10 2 km 3 for the model and CARS, respectively; Fig. 13).Contrary to the hypoxic volume, the suboxic volume is affected by parameter values associated with the nitrification process.Indeed, the suboxic volume is divided by 1.6 when K N42 and K N23 are divided by 10 (see the volume change between Tests no. 4 and 5 and between Tests no. 14 and 15 in Fig. 13).In the OMZ, NO − 3 reduction generates a NO − 2 source that is quickly oxidized to NO − 3 by the second stage of nitrification.This link between both processes in the OMZ generates a maximum of nitrifying activity (which uses O 2 to oxidize NH + 4 and NO − 2 to NO − 3 ) directly affecting the suboxic volume.However, it has no impact on the hypoxic area as NO − 3 reduction process is not present for oxygen concentrations above 25 mmol O 2 m −3 .So the tested parameter values associated with organic matter decomposition, vertical sinking, and nitrification directly control the OMZ volume, but interestingly have no influence on averaged and minimum oxygen concentrations in the OMZ.
The "reference simulation" gives satisfying results as compared to CARS oxygen volume, except on the very coastal area.None of our simulations can increase the simulated OMZ volume close to the coast (bathymetry shallower than 130 m).Coastal sediment processes (a supplementary oxygen demand) play a key role in maintaining bottom oxygen concentrations at low values over the continental shelf, and they will be implemented in a future work.However, impacts of sediment biogeochemistry on water column N-C biogeochemical cycles are still uncertain.Bianucci et al. (2012) studied the role of sediment denitrification in water column oxygen dynamics of the northern end of the California Current system, and argued that sediment denitrification is balanced by water-column nitrification resulting in bottom oxygen concentrations equivalent to the ones simulated without sediment module.Conversely, in regions where recycled production is important, losses of bioavailable nitrogen through sediment denitrification can decrease the primary production and then increase the bottom oxygen concentrations (Bianucci et al., 2012).
The simulated OMZ off Namibia does not extend more poleward than 26 • S (Fig. 14a), and represents the southern boundary of the eastern tropical southeast Atlantic zone (Monteiro and van der Plas, 2006), supplied by the poleward undercurrent on the shelf and slope of the Benguela.Simulated minimum oxygen concentrations are around 20 mmol O 2 m −3 on average over the 8 analyzed years (Fig. 14a), but lower concentrations can be simulated on 3-day averages (not shown; its occurrence can be seen in Fig. 9c).The minimum oxygen concentrations are located on the bottom waters over the continental shelf, and remain around 250-300 m depth off the shelf break (Fig. 14b).On 8 yr average, the OMZ thickness is about 20 m over the shelf break (Fig. 14c), but it can increase up to 150 m over the slope when the simulated water column is particularly deoxygenated (using 3-day averages; not shown).Monteiro and van der Plas (2006) showed that the remote and local forcing sparked off low oxygen variability in the Namibian subsystem.The input fluxes at northern and southern boundary conditions govern the OMZ formation, while local biogeochemical oxygen demand only acts on persistence and intensity of the OMZ (Monteiro and van der Plas, 2006).The 3-D coupled modeling experiment presented here allows an analysis of the controlling processes of this OMZ.The 8 analyzed years present a balanced oxygen budget.Physical and/or biogeochemical processes only act to conserve these low oxygen concentrations.Our analysis confirms that oxygen advection is the main sink term maintaining the OMZ offshore the isobath 300 m whereas mixing and mainly its vertical contribution act to dissipate the OMZ (oxygen source; Fig. 14d and e).This finding is also consistent with a modeling study covering the whole Benguela system using the SAfE configuration with the same coupled ROMS/BioEBUS model (Le Vu et al., 2013).It also highlights that the minimum oxygen concentrations located between isobaths of 100 and 300 m are not driven by physical forcing (Fig. 14d) as in this area; advection processes generate a weak gain of oxygen, when considering the 8 yr average.Instead, the spatially constant biogeochemical sink (Fig. 14f), mainly due to nitrifying bacteria, becomes the only process to deplete the oxygen content between isobaths 100 and 300 m.This spatial shift in governing processes maintaining the minimum oxygen concentrations occurs over the whole analyzed domain, from 20 • S to the southern limit of the OMZ off Namibia, and has not been reported previously to our knowledge, and has to be confirmed with future studies.Supplementary sensitivity tests (not shown) using different O 2 concentrations at the northern boundary (here 19 • S) were also performed.These tests confirmed that meridional advection of oxygen-depleted waters at this northern boundary represents a necessary condition to form the OMZ off Namibia but does not drive its intensity, confirming the results of Monteiro and van der Plas (2006).
In this work we investigated the importance of the biogeochemical parameter values on the OMZ features, using an already validated configuration for the physics in the studied area and the associated physical settings (Veitch et al., 2009).However, previous works have shown that physical parameter values such as the mixing coefficients can also have a great impact on the OMZ characteristics (i.e., Duteil and Oschlies, 2011;Banyte et al., 2012;Fischer et al., 2012).
OMZs associated to EBUS also play a major role in the global nitrogen cycle acting as areas of N losses as well as N 2 O production and outgassing.We now discuss the impact of parameter values on these two mechanisms.

N losses
Denitrification and anammox processes reduce bioavailable nitrogen to gaseous end products.These losses of fixed inorganic nitrogen from the marine ecosystem were estimated using the coupled model, and analyzed through the sensitivity analysis on biogeochemical parameter values (Fig. 13).
Due to its low concentrations NO − 2 is here the limiting factor for both processes (denitrification and anammox).As a consequence, increasing denitrification rates (Test no. 9 in Fig. 13) directly impact N losses as NO − 3 reduction supplies the NO − 2 pool for the second stage of denitrification (NO − 2 reduction) that generates N loss.On the other hand, increasing anammox rate alone cannot initiate significant N losses.To increase N losses due to anammox, processes supplying the NO − 2 pool (NH + 4 oxidation and NO − 3 reduction) have to be more efficient.Thus, slower nitrifying activity associated with increased anammox and denitrification rates (Test no.14) improves the anammox fluxes off Namibia (see Sect. 4.2).At the same time, it increases N losses due to the anammox process more than N losses due to denitrification (Fig. 13), in agreement with Kuypers et al. (2005) off Namibia and Thamdrup et al. (2006) off northern Chile.However, this combined set of parameters also affects the N 2 O production in the water column off Namibia (see Sect. 4.3).We could also independently change the nitrification rates (accelerate the NH + 4 oxidation rate and slow down the NO − 2 oxidation rate) in order to enhance NO − 2 pool available for anammox bacteria.However, daily rates are close to in situ data (Fig. 11), and furthermore higher NO − 2 oxidation rates than NH + 4 oxidation rates are observed in the data due to significant source of nitrite by NO − 3 reduction (Lam et al., 2009;Füssel et al., 2011).In Test no. 15 we increase the anammox and denitrification rates without changing the nitrifying activity.This combination results in NO − 2 reduction as the main process to deplete the water masses of the Namibian upwelling system in nitrogen.
On the 8 yr mean, the anammox process only contributes to one fifth of the total N losses (Fig. 13).This finding is not in agreement with Kuypers et al. (2005) and Thamdrup et al. (2006), who found that the dominant N 2 production mechanism was anammox in the OMZ off Namibia and northern Chile, respectively, using punctual observations during the same autumn season (Namibia: RV Meteor cruise in March/April 2003 -4 stations; Chile: DINAMO cruise in March 2004 -2 stations).However, anaerobic ammonium oxidation to N 2 using nitrite as an oxidant does not seem to depend on other microorganisms for the production of nitrite and ammonium as it is currently parameterized in biogeochemical models.Instead anammox bacteria can reduce nitrate to ammonium via dissimilatory nitrate and nitrite reduction (Kartal et al., 2007).The apparent complexity of anammox bacteria deserves revision of their representation and formulation in the models, thereby increasing in situ measurements.
It has become clear that anammox bacteria are of global importance in the biogeochemical nitrogen cycle (Strous and Jetten, 2004;Arrigo, 2005).It is currently assumed that anammox process contributes up to 50 % to the global fixed nitrogen removal (Arrigo, 2005); nevertheless, the relative importance of anammox and denitrification processes in N losses is still an open discussion due to contrasting results in OMZs.

N 2 O concentrations and emissions to the atmosphere
Using our sensitivity analysis on N 2 O parameterizations, we choose the one from Suntharalingam et al. (2000Suntharalingam et al. ( , 2012)), with a higher β coefficient (see Sect. 4.3), and implemented it in the "reference simulation".Indeed, between the "standard simulation" using the parameterization of Nevison et al. (2003) (Test no. 4) and the "reference simulation" using the one from Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) (Test no.15), mean N 2 O concentrations increased by a factor of 20 % in the intermediate layer (100 to 400 m depth) over the continental shelf and shelf break where the OMZ is developed (Fig. 12a).This increase provides improvement in maximum N 2 O concentrations within the range of the first observations available in the Namibian upwelling system (see Sect. 3.5 and Fig. 12b).The importance of parameter values and N 2 O parameterization becomes even evident when zooming in on the OMZ area as mean N 2 O concentrations increased by a factor of 50 % between the "standard" and the "reference simulations" (Fig. 13).The OMZ contributes to 40.5 % of total N 2 O production over the analyzed domain (22-24 • S and 10-15 • E) and 25 % over the whole simulated domain, in agreement with Suntharalingam et al. (2000).Indeed, these authors suggest that 25 to 50 % of the global oceanic N 2 O source could be localized in low-oxygen regions.
The parameter values and N 2 O parameterization directly impact N 2 O concentrations in the OMZ, but do not significantly constrain the N 2 O sea-air fluxes.Indeed, N 2 O outgassing only increases by 5.5 % between both simulations (standard and the reference simulations; Fig. 13) due to similar mean N 2 O concentrations in the mixed layer (increase < 1 %; Fig. 13).These results diverge from those of Suntharalingam et al. (2000).They found that the majority of N 2 O formed under low-oxygen conditions is directly outgassed to the atmosphere, suggesting that an increase of N 2 O production in the OMZ should result in similar increase of N 2 O sea-air flux.In our simulation, N 2 O production associated with the OMZ of the Walvis Bay area is mostly advected southward by the poleward undercurrent instead of being locally outgased to the atmosphere (Gutknecht et al., 2013).Our budgets especially underestimate the N 2 O production inshore of 130-m isobath as the OMZ is not correctly reproduced on the continental shelf.
Improvement of parameter values associated with organic matter decomposition, vertical sinking, and nitrification appears as important as the parameterization used to estimate N 2 O emissions.Indeed, between the "standard simulation" using the Nevison et al. (2003) parameterization (Test 4) and the "reference simulation" using the one from Suntharalingam et al. ( 2012) with β = 0.03 mol N 2 O mol N −1 (Test 15), 45 % of N 2 O emission increase comes from the parameterization used in BioEBUS and about 55 % due to the improvement of the parameter values.So, both improvements are crucial to better simulate N 2 O emissions using oceanic biogeochemical models.However, in situ measurements are needed to refine the α and β coefficients used in the parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) and to improve the N 2 O vertical distribution and emissions to the atmosphere.
Recently, Zamora et al. (2012) suggested that N 2 O production increases linearly with decreasing oxygen, rather than exponentially, in the eastern tropical Pacific, based on MEMENTO database.This fundamental difference in N 2 O production parameterization introduces a large uncertainty in biogeochemical models, especially in EBUS and associated OMZs (Zamora et al., 2012), even more in a context of expansion and intensification of the OMZs (Stramma et al., 2008;Deutsch et al., 2010).Moreover, it is currently assumed that N 2 O is produced via two mechanisms: denitrification under suboxic conditions with N 2 O as an obligatory gaseous intermediate to N 2 , and nitrification in oxic environments with N 2 O as a by-product.However, NO detoxification by anammox could be another source of N 2 O (Strous et al., 2006;Kartal et al., 2007).So more studies, especially derived proxies for existing communities (e.g., ladderanes as in Kuypers et al., 2005) are needed to better understand the N 2 O cycle and improve its representation in biogeochemical models, however our results are promising.
Finally, we currently do not explicitly account for the consumption of N 2 O during the second step of denitrifica-tion as the parameterization of N 2 O production of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) is only based on O 2 concentrations.N 2 O consumption occurs at very low O 2 concentrations (< 1-2 mmol O 2 m −3 ; Gruber, 2004) or maybe up to 10 mmol O 2 m −3 , as reported by Zamora et al. (2012) in the eastern tropical Pacific.However, N 2 O consumption is not a relevant process in our "reference simulation"; oxygen minimum values usually remained above the threshold limit of 10 mmol O 2 m −3 (see 3-day averaged oxygen concentrations in Fig. 9c).We will include this consumption process using the formulation of Jin and Gruber (2003) or Zamora et al. (2012) in future work.

Conclusion
We developed a biogeochemical model (BioEBUS) taking into account the main processes linked with EBUS and associated OMZs.We implemented a 3-D realistic coupled physical/biogeochemical configuration for the Namibian upwelling system using the high-resolution hydrodynamic ROMS model.A detailed validation using in situ and satellite data as well as "metrics" is presented in this paper.Sensitivity analyses were performed on key parameters and different N 2 O parameterizations to adjust the coupled model, and relevant simulations are discussed in order to improve the spatio-temporal distribution of oxygen, nitrate, Chl a, and N 2 O concentrations, as well as key processes of the microbial loop.The "reference simulation" gives satisfying results as compared to in situ and satellite data as well as microbial rates available in the studied area.We then discussed the impact of parameter values on the OMZ off Namibia, N losses due to denitrification and anammox processes, and N 2 O concentrations and emissions to the atmosphere.In a companion paper (Gutknecht et al., 2013), this same configuration is used to investigate the full N budget off Walvis Bay in order to understand and quantify the N offshore export likely to sustain the primary production of the subtropical gyre, the export production and subsequent losses of fixed N via denitrification and anammox under suboxic conditions, and the N 2 O outgassing in the upwelling area.
The low-oxygen water masses are very sensitive to biogeochemical parameter values associated with organic matter decomposition, vertical sinking, and nitrification.After adjustment of the coupled model, the simulated hypoxic water volume is similar to the observations (CARS climatology, 2006), and the suboxic volume also compares well with observations when excluding the area shoreward of the 130 m isobath.Coastal sediments appear necessary to improve the oxygen distribution in the coupled model.
Using the high-resolution configuration, a spatial decorrelation of dominant processes maintaining the minimum oxygen concentrations off Namibia was highlighted with the poleward meridional advection of oxygen-depleted waters offshore of the 300 m isobath and nitrifying activity inshore of this isobath.
After sensitivity analyses that clearly allowed improving NO − 3 reduction and anammox fluxes and associated N losses, NO − 2 reduction is the main process to deplete the water column of the Namibian upwelling system in nitrogen; anammox process only contributes to 20 % of the total N losses.This relative contribution is lower than the currently assumed contribution (up to 50 %) for the global ocean.However, this relative contribution of anammox and denitrification processes in oceanic N losses still needs to be better understood as anammox bacteria appear to have a more complex and versatile behavior than previously believed.Some biogeochemical pathways are not yet well understood in natural ecosystems.The sensitivity analysis also highlights the importance of nitrification processes and the need to take explicitly into account NO − 2 in such low-oxygen environments.Both steps of nitrification -ammonium oxidation and nitrite oxidationare necessary to improve the representation of microbial activity linked with the OMZ, and should be detailed in coupled models dedicated to EBUS.
The N 2 O parameterization of Suntharalingam et al. (2000Suntharalingam et al. ( , 2012) ) using a higher β coefficient gives satisfying results as compared to in situ measurements made off Namibia, and especially for mean N 2 O concentrations in the OMZ.The OMZ area off Namibia contributes to 25-40 % of total N 2 O production in the studied area; the same relative contribution as for global estimations.However, N 2 O concentrations in the OMZ do not directly constrain the N 2 O sea-air fluxes as N 2 O production associated with the OMZ is mostly advected southward by the poleward undercurrent.
Lastly, a large number of uncertainties still remains on the understanding of biogeochemical cycles, and especially the nitrogen cycle.For instance, the role of the anammox process in dissimilatory nitrate reduction is still to be elucidated.N 2 O production parameterization introduces a large uncertainty in biogeochemical models, especially in EBUS and associated OMZs (Zamora et al., 2012).These questions are of first importance to better understand the nitrogen cycle, and better apprehend the expansion and intensification of the OMZs already observed in a context of climate change.Therefore, more studies are needed to better understand the N cycle and improve its representation in biogeochemical models.

Fig. 2 .
Fig. 2.Interactions between the different compartments of the BioEBUS model.Black arrows represent the nitrogen-dependent processes; red arrows, the oxygen-dependent processes; blue arrows, the processes linked with N 2 O production.To simplify the representation of all interactions between variables, arrows from or to a grey rectangle act on all variables included in this grey rectangle.For example, the arrow between nutrients and phytoplankton (assimilation) is a simplification of 6 interactions: NO − 3 to P S , NO − 3 et al.: Coupled physical/biogeochemical modeling including O 2 -dependent processes microphytoplankton concentrations, dz is the depth step (vertical thickness; in m).

Fig. 3 .
Fig. 3. Cruise sections and stations used for the model/data comparison: the Danish Galathea expedition (October 2006), the ME-TEOR expedition 57/2 (February 2003), the AHAB1 expedition (January 2004), AMT 6 cruise (May 1998), the Walvis Bay mooring at 23 • S between 1994 and 2004 (Monteiro and van der Plas, 2006), the time series of copepod abundance from the coast (14.5 • E) to 70 nautical miles (13.23 • E) off Walvis Bay at 23 • S between 2000 to 2007 (Kreiner and Ayon, 2008), and the N 2 O data from FRS Africana cruise (December 2009).The black dashed lines indicate the bathymetry (in meters).Specific stations are indicated and will be discussed later in the text.

Fig. 4 .
Fig. 4. Taylor diagrams for (a) temperature, (b) salinity, (c) density, (d) oxygen, and (e) nitrate concentrations.The radial distance from the origin is proportional to the standard deviation of a pattern (normalized by the standard deviation of the data: in situ, satellite or climatological data).The green dashed lines measure the distance from the reference point (the data; black filled circle in the diagram) and indicate the RMS error.The correlation between both fields is given by the azimuthal position.The statistics use data from the METEOR expedition 57/2 in February 2003 (transects at 23 • S, 24.4 • S, and 25.5 • S), the AHAB1 expedition in January 2004 (transects at 22 • S, 23 • S, 24 • S, 25 • S, 26 • S, and 27 • S), the Galathea data in October 2006 (surface data, 5 stations, and triaxus data), the AMT 6 cruise in May 1998 (transect at 26.7 • S and 4 stations), and the CARS databases(2006, 2009)  (seasonal and annual climatologies).For the comparison with in situ data, simulated fields were monthly averaged (note there are 10 outputs by month) using the last 8 yr of simulation (Y12-Y19).For the comparison with CARS database, simulated fields were either seasonally or annually averaged.

Fig. 5 .
Fig. 5. Simulated 8 yr mean of (a) temperature ( • C), (b) salinity, (c) density (kg m −3 ), (d) oxygen (mmol O 2 m −3 ) and, (e) nitrate (mmol N m −3 ) concentrations at 23 • S between the surface and 600 m depth, with a zoom in the first 100 m depth.Colored circles for the annual mean CARS databases (2006, 2009) are overlaid on the simulated fields using the same color bar.

Fig. 7 .
Fig. 7. 8 yr mean and times series (last 8 yr of simulation averaged between 22-24 • S) of logarithm of surface chlorophyll a concentrations (mg Chl m −3 ) from the coupled model (a, c) and SeaWiFS climatology (b, d).Taylor diagram for logarithm of chlorophyll a concentrations (e).See legend of Fig. 4 for details.Biases between the modeled fields and compared data have been added to the right of the legend.

Fig. 8 .
Fig. 8. Vertical distribution of Chl a(mg Chl m −3 ) at stations 14 (a), 15 (b), and 16 (c).The black solid line represents the simulated field (averaged May over 8 yr).The markers represent the estimations made during the AMT 6 cruise, in May 1998, from four different sources (E.Fernandez for the diamonds, P. Holligan for the circles, the British Oceanographic Data Centre for the stars, and R. Barlow for the triangles).

Fig. 9 .
Fig. 9. (a) Oxygen (mmol O 2 m −3 ), and (b) nitrous oxide (10 −3 mmol N 2 O m −3 ) concentrations estimated with the coupled model at 23 • S (averaged December over 8 yr) between the surface and 1500 m depth, with a zoom in the first 100 m depth.Colored circles for the FRS Africana (December 2009) data are on the modeled fields using the same color bar.(c) Nitrous oxide concentrations (10 −3 mmol N 2 O m −3 ) as a function of oxygen concentrations (mmol O 2 m −3 ) at 23 • S. Black dots represent the simulated 3-day averaged fields for the whole analyzed period (Y12-Y19), and red stars in situ measurements from the FRS Africana cruise in December 2009.

Fig. 10 .
Fig. 10.Eight-year averaged vertical profiles for O 2 (a) and NO − 3 (b) concentrations from the surface to 600 m depth, and for chlorophyll a concentrations (c) from the surface to 100 m depth, for different tests.Vertical profiles were averaged between 22-24 • S and 10-15 • E. For oxygen, vertical profiles excluding the continental shelf in (a).
Fig. 11.NH + 4 and NO − 2 oxidation rates (10 −2 mmol N m −3 d −1 ) in the whole water column (top) and in the OMZ (middle).Bottom: NO − 3 reduction and anammox rates (10 −2 mmol N m −3 d −1 ).Rectangles represent 8 yr mean values for the different tests performed for the sensitivity analyses, averaged between 22-24 • S and 10-15 • E, and black bars represent the standard deviation of the 8 annual mean values (Y12-Y19) for these different quantities.The range of in situ data is added to the right of the histogram (values are reported in Table6).

60Fig. 13 .
Fig. 13.Top: oxygen volume below 25 mmol O 2 m −3 (10 2 km 3 ; suboxic volume or OMZ volume) and below 60 mmol O 2 m −3 (10 3 km 3 ; hypoxic volume).The same fields using CARS (2006) climatology are added to the right of the histogram.Middle: N loss due to denitrification and anammox processes (10 8 mol N yr −1 ).Bottom: mean N 2 O concentrations (mmol N 2 O m −3 ; left axis) in the mixed layer (∼ 0-50 m) and in the OMZ, and mean N 2 O emissions to the atmosphere (10 8 mol N yr −1 ; right axis).There are 8 yr mean values for the different tests realized for the sensitivity analyses, computed between 22-24 • S and 10-15 • E. Black (or blue for sea-air flux) bars represent the standard deviation of the 8 annual mean values (Y12-Y19) for these different quantities.

Fig. 14 .
Fig. 14.(a) Minimum oxygen concentrations (mmol O 2 m −3 ) simulated in the OMZ, (b) depth (m) of these minimum oxygen concentrations, (c) thickness (m) of the OMZ, (d) oxygen advection term (sum of zonal, meridional and vertical components; mmol O 2 m −3 d −1 ;), (e) oxygen mixing term (sum of horizontal and vertical components; mmol O 2 m −3 d −1 ), and (f) oxygen consumption by biogeochemistry (oxic decomposition + nitrification; mmol O 2 m −3 d −1 ) at the depth of minimum oxygen concentration.For fluxes, negative values represent an oxygen sink and positive a source.All fields are averaged over the 8 analyzed years.

Table 1 .
BioEBUS biogeochemical model state variables (symbols and units), initial surface values used for the simulation, and scale depth considered for the exponential decrease with depth (for NO − 2 , NH + 4 and DON).
than in situ data and CARS climatology (2009) at 23 • S (not shown) and Chl a concentrations from WOA (2001).

Table 5 .
Sensitivity analyses: tested parameter values, N 2 O parameterizations, and combined simulations.