Tidal impacts on primary production in the North Sea

This study highlights the importance of tides in controlling the spatial and temporal distributions of phytoplankton 5 and other factors related to growth, such as nutrients and light availability. To quantify the responses of net primary production (NPP) to tidal forcing, we conducted scenario model simulations considering M2 and S2 tidal constituents using the physicalbiogeochemical coupled model ECOSMO. The results were analysed with respect to a reference simulation without tidal forcing, with particular focus on the spatial scale of the tidally induced changes. Tidal forcing regulates the mixingstratification processes in shelf seas such as the North Sea and hence also influences ecosystem dynamics. In principle, the 10 results suggest three different response types with respect to primary production: i) in southern shallow areas with strong tidal energy dissipation, tidal mixing dilutes phytoplankton concentrations in the upper water layers and thereby decreases NPP. Additionally, tides increase turbidity in near-coastal shallow areas, which has the potential to further hamper NPP. ii) In the frontal region of the southern North Sea, which is a transition zone between stratified and mixed areas, tidal mixing infuses nutrients into the surface mixed layer and resolves summer nutrient depletion, thus sustaining the NPP during the summer 15 season after spring bloom nutrient depletion. iii) In the northern North Sea, the NPP response to tidal forcing is limited. Additionally, our simulations indicate that spring bloom phenology is impacted by tidal forcing, leading to a later onset of the spring bloom in large parts of the North Sea and to generally higher spring bloom peak phytoplankton biomasses. By testing the related changes in stratification, light conditions and grazing pressure, we found that all three factors potentially contribute to the change in spring bloom phenology with clear local differences. Finally, we also analysed the impact of the spring-neap 20 tidal cycle on NPP. The annual mean impact of spring-neap tidal forcing on NPP is limited. However, locally, we found substantial differences in NPP either in-phase or anti-phase with the spring-neap tidal cycle. These differences could be attributed to locally different dominant factors such as light or nutrient availability during spring tides. In general, we conclude that in shallow shelf seas such as the North Sea, intensified vertical mixing induced by tidal forcing could either promote NPP by counteracting nutrient depletion or hinder NPP by deteriorating the light environment because of the resuspension and 25 mixing of suspended matter into the euphotic zone.

focused on the North Sea region between 5°W-9.5°E and 48°N-58.5°N because tides are only of minor relevance in the Baltic Sea. To resolve thermal stratification in the upper water column, the vertical resolution was set to 5 m in the upper 40 m of the 95 water column and decreased gradually with depth below 40 m. To reduce numerical diffusion in the implemented upwind advection scheme, a shape-preserving total variation diminishing (TVD) scheme (Yee et al., 1985) was adopted, which significantly improved the representation of hydrodynamics and ecosystem processes, especially processes related to fronts.
A detailed description of the method and model responses to the changed advection scheme has been provided by Barthel et al. (2012). 100 The biogeochemical component of ECOSMO was developed to describe the lower trophic level dynamics of the marine ecosystem using a N(utrient)P(hytoplankton)Z(ooplankton)D(etritus) conceptual model framework. The ecosystem model component was first introduced for the North Sea by Schrum et al. (2006a) and further developed for a wider range of ecosystems, including relevant characteristics for the Baltic Sea, by Daewel and Schrum (2013). Detailed validations against nutrient observations have shown that the model is capable of simulating lower trophic level ecosystem dynamics in the North 105 Sea, and the temporal variability at inter-annual to decadal scales simulated by ECOSMO could be corroborated by observations (Daewel and Schrum, 2013). ECOSMO simulates the nutrient cycling of silicate, phosphorus and nitrogen in the water column and in the sediments considering processes such as primary production, grazing and excretion by zooplankton, remineralization and sediment-water coupling. A detailed description of the ecosystem model is given in Daewel and Schrum (2013). In total, 16 state variables were solved, including three functional groups for primary producers (diatoms, flagellates 110 and cyanobacteria). In the second trophic level, two groups of zooplankton were considered and differentiated based on feeding preferences. To additionally account for the shading effects of DOM and detritus, which were not considered in Daewel and Schrum (2013), the formulation of light attenuation was modified as previously suggested by Nissen (2014). To capture the productive and turbid characteristics, dissolved organic matter (DOM) was parameterized by fast remineralization rates and a low sinking velocity, in contrast to the fast sinking velocity and slow remineralization rates of particulate organic matter 115 (detritus). Therefore, the vertical light attenuation consisted of background attenuation ( 1 ) (induced by the water body and inorganic SPM), phytoplankton self-shading ( ), and additional shading impacts of DOM ( ) and detritus ( ), as shown in Eq.1.
(1) 120 While background attenuation 1 (0.03 −1 , Urtizberea et al., 2013) remained constant in the water column, self-shading depended on both (0.2 2 −1 ) and the phytoplankton concentration (P). As suggested by Stedmon et al. (2000) and Tian et al. (2009), and detritus were set to 0.29 2 −1 and 0.2 2 −1 , respectively. Compared to Daewel and Schrum (2013), these changes enabled the dynamical coupling of turbidity to the seasonal production cycle, as previously 125 discussed by Nissen (2014). A corresponding validation of surface nutrients and comparison of mean primary production (Appendix B) confirms that the performance of ECOSMO in the North Sea region changes only marginally with respect to the original model version. Frontal production and production in deeper stable stratified waters increased slightly, while production near the coast was slightly decreased. The model is thereby capable of resolving tidal influences on primary production via potentially competing processes. Tidal mixing releases nutrient limitation, thus fostering NPP, but tides also cause the 130 resuspension and mixing of suspended matter into the euphotic zone, which reduces light availability in the water column, thus reducing NPP. In addition to relevant bottom-up processes, the model also resolved phytoplankton-zooplankton feedbacks and vertical oxygen and temperature profiles, which alter the remineralization of organic matter and consequently nutrient cycling.
Besides organic matter contribution to light shading, in the coastal area inorganic SPM also has the potential to filter light and reduce primary production. We do not consider a dynamic coupled SPM modelling approach, but consider a simplified 135 consideration of inorganic SPM through implementing the background attenuation. To address the uncertainties related to SPM, we tested the effect of inorganic SPM on our findings with help of an additional numerical simulation, where we implemented a climatological SPM field (daily resolution, with 31 vertical layers in the original dataset) (Große et al., 2016;Heath et al., 2002) and added the SPM's contribution to the light attenuation scheme. Details of the inorganic SPM data set and implementation are given in the Appendix C. The results confirmed the validity of our assumption that the spatial 140 variability of SPM can be neglected for the sensitivity study performed here. Despite the existing effect of inorganic SPM on light conditions and spatial variability of inorganic SPM, there is only minor sensitivity found for the case studies tidal vs nontidal forcing, and Eq.1 can be considered as a proper parameterisation within the context of our study.
The ability to properly resolve intensified frontal production and the consideration of key processes influencing light and nutrient limitation related to tidal forcing make ECOSMO an appropriate tool to assess tidal impacts on NPP in the spatially 145 highly diverse North Sea. As already stated in Daewel and Schrum (2013), ECOSMO estimates of annual NPP in the North Sea (Fig.2) are at the lower edge of what has been simulated for the area (Holt et al., 2012;van Leeuwen et al., 2013). The relatively low estimates mainly appear in the northern North Sea (NNS), where primary production is estimated to be approximately 125 −2 −1 based on observations (Van Beusekom and Diel-Christiansen, 1994), and in the European continental coast, where NPP observations range between 199-261 −2 −1 (Joint and Pomroy, 1993). The simulation fits 150 well with observation-based estimates of NPP in the British coast of approximately 75-79 −2 −1 (Joint and Pomroy, 1992) and primary production estimates of 100 −2 −1 and 119-147 −2 −1 in the central parts of the North Sea and at Dogger Bank, respectively (Joint and Pomroy, 1993).

Model setup
A detailed description of the model setup was given by Daewel and Schrum (2013); therefore, we will only provide a brief 155 overview of the forcing data used for the model simulation, particularly emphasizing the changes made to the previously described setup. These changes mainly concern the river discharge and nutrient load data sources. The simulation was initialized in 1948 using climatological data from the World Ocean Atlas (WOA) (Conkright et al., 2002) for nutrients and observational climatology for temperature and salinity (Janssen et al., 1999). The full simulation period encompasses 68 years, ending in 2015, and is forced with atmospheric boundary conditions provided by the NCEP/NCAR reanalysis (Kalnay et al., 160 1996). Additional forcing data include wet deposition for nitrogen, which were prescribed using data from a Community Multiscale Air Quality (CMAQ) model (Matthias et al., 2008), and boundary values for nutrients, temperature and salinity at the open boundaries to the North Atlantic, for which we used the same climatological data as those used for the initial conditions. For salinity, additional annual anomalies were retrieved from observational data available at the ICES (International Council for the Exploration of the Sea) database (http://www.ices.dk). An updated set of river runoff and nutrient 165 load data was applied with more complete river forcing data coverage for the North Sea and Baltic Sea. A multitude of data were provided by Sonja van Leeuwen (Royal Netherlands Institute for Sea Research, pers. communication) containing the following datasets: UK data were processed from raw data from the Environment Agency (England and Wales, contains Natural Resources Wales information © Natural Resources Wales and database rights), the Scottish Environment Protection Agency (Scotland), the Rivers Agency (Northern Ireland) and the National River Flow Archive. French water quality data 170 were provided by Agence de l'eau Loire-Bretagne, Agence de l'eau Seine-Normandie, OSUR web Loire-Bretagne and SIEAG (Systeme d'information sur l'eau du bassin Adour Garonne), while daily flow data were obtained from Le Banque Hydro (www.hydro.eaufrance.fr). German and Dutch riverine data were provided by the University of Hamburg (Pätsch, J., Lenhart, 2004). Norwegian water quality data were provided by the Norwegian Water Resources and Energy Directorate (NVE), with daily flow data supplied by the Norwegian Institute for Water Research (NIVA). Danish water quality data were provided by 175 the National Environmental Research Institute (NERI). Water quality data for Baltic rivers were provided by the University of Stockholm and the Baltic Nest. Furthermore, nutrient status and freshwater runoff information in the southern and eastern Baltic Sea was supplemented by data from the Balt-Hype model (Arheimer et al., 2012;Lindström et al., 2010). Nutrient loads from Danish waters were provided by Marie Maar (pers. communication) and were similar to the forcing data used for the HBM-ERGOM simulation (Maar et al., 2016). These data stem from a national monitoring program (Windolf et al., 2011(Windolf et al., , 180 2012) and from the hydrological Denmark model, which provides runoff calculations for ungauged areas of Denmark (Henriksen et al., 2003).
We selected a relatively short time period  for our analysis to assure a long enough spin-up time that accounts for the characteristic long time scales of the North Sea-Baltic Sea system (Daewel and Schrum, 2013). The period from 1990-2015 will hereafter be called the analysed period. Tidal cycles with long periods, such as the nodal and elliptical cycles, 185 although considered in the forcing via nodal corrections of partial tide amplitudes and phases (see 2.3), are not targeted in this study.

Tidal forcing and scenarios
Sea surface elevation was prescribed at the open boundaries, with a time step of 20 min. Daily mean sea surface elevation data were taken from a diagnostic model simulation for the wider northeast European Shelf (Backhaus and Hainbucher, 1987) and 190 also forced with the NCEP/NCAR reanalysis. In addition, tidal elevations were calculated from tidal constituents provided by the German Federal Maritime and Hydrographic Agency (Federal Maritime and Hydrographic Agency Deutsches Hydrographisches Insititut, 1967) . Nodal corrections were implemented in the calculation of tides to represent the long-term variation in lunar nodes. For the standard tidal scenario, partial M2 tide (principle lunar tide) and S2 tide (principle solar tide) (Thomson and Emery, 2014) were considered; we hereafter call this scenario the tidal scenario. To evaluate the contribution 195 of the spring-neap tidal cycle, a tidal scenario using only the M2 partial tide, called the M2 scenario, was simulated and discussed in comparison to the tidal scenario. To quantify the overall impact of tidal forcing, a scenario without tidal forcing at the open boundary was simulated to yield the non-tidal reference state of the system (non-tidal scenario).

Postprocessing of model results
The responses of ecosystem productivity to tidal forcing were assessed by comparing the annual mean NPP during the analysed 200 period between the tidal and non-tidal scenarios (tidal scenario minus non-tidal scenario). Furthermore, we disentangled processes that might contribute to variations in NPP, such as the seasonality of spatial patterns in limitation factors (nutrients vs. light), spring bloom phenology, the impacts of the spring-neap cycle on NPP variability, and the contribution of subsurface production to the overall NPP. We quantified these processes using subdomains and further made comparisons between scenarios, emphasizing spatial variability and the seasonal cycle. 205

Subdomain division and identification of representative grid cells for process based analysis.
The pre-division of the area into subdomains is based on a combination of geographic location, bathymetry and the local responses of NPP to tidal forcing (increase, decrease). First, SNS and NNS were divided by the 65 m isobath. In the SNS, areas with positive and negative NPP response to tides were separated (Fig.2). The negatively responding area in the SNS was further geographically divided into the English Channel (EC, south of 52°N) and an area along the continental coast (neg. 210 SNS). In the NNS, the area of the Norwegian Trench (NT) characterized by a water depth deeper than 200 m was separated.
The remaining region of the NNS was further divided based on the response of NPP to tidal forcing. The area along the eastern British coast (BC), which shows elevated NPP in response to tides was separated from the negative responding area in the middle of NNS (deep NNS) (Fig.2). In the east of the NNS, a separate area with mild increase of NPP was identified (low-sen. NNS). Based on this pre-division of subdomains, we identified the most representative grid cell within each subdomain using 215 correlation analysis (Eliasen et al. 2017,) (Fig.A1). To identify the most representative grid cell location in each subdomain, we first produced a time series of the NPP differences between the non-tidal scenario and tidal scenario for each grid cell.
Subsequently, we estimated, for each of the grid cells, the correlation to the time series of the other grid cells within the same pre-divided subdomain. The grid cell with the highest correlation coefficient to all other grid cells in each subdomain was selected as the most representative point for further analysis. 220

Quantification of key processes controlling the spring bloom
The peak amplitude and the onset time of the spring bloom for the different scenarios were compared. The onset of the spring bloom is defined here as the day when the daily vertically integrated NPP reaches its maximum prior to the spring maximum in diatom biomass (Fig.A2) (Sharples et al., 2006). Diatom time series were pre-processed by a 15-day running mean to remove short-term maxima induced by the spring-neap tidal cycle (Sharples et al., 2006). To further disentangle mechanisms resulting 225 in spring bloom phenology differences among the scenarios, we quantified potentially related biological and physical factors relevant for spring bloom dynamics, such as the zooplankton biomass prior to the onset of the spring bloom, light conditions and development of stratification, for each grid cell.
In particular i) the vertically averaged zooplankton biomass in the winter season (January & February) was considered a proxy for potential grazing pressure at the beginning of the growth season. ii) The integrated value of the light-limiting term in the 230 upper 50 metres of the water column was used to estimate the light conditions for phytoplankton growth. To quantify the time when the light was sufficient for phytoplankton growth in each year, we estimated the date when the integrated light limiting term exceeded 0.85 for three consecutive days. iii) The stratification was recognized as critical temperature difference (ΔT) between surface layers and layers below exceeding 0.5 ℃ . Similar methods have been used in many other studies (Gong et al., 2014;Karl and Lukas, 1996;Richardson et al., 2002). For the identification of the onset time of stratification, ΔT exceeding 235 0.5°C for three consecutive days is required. The time window (3 days) (Sharples et al., 2006) was chosen to filter out shortlived stratification variations and the day-night heating/cooling cycle. The mixed layer depth is defined as the thickness of the surface mixed layer, ranging from surface to pycnocine. iv) The averaged mixed layer depth in May was used as a measure for stratification depth.
The onset of the spring bloom, the first day of the year with stratification and the first day of the year with sufficient light 240 conditions were identified for each grid point for every simulated year; subsequently the percentage of years in which those time identifiers were advanced or delayed in the tidal scenario compared to that in the non-tidal scenario as a response to tidal forcing was estimated for every grid cell. The tidal induced increase/decrease of winter zooplankton biomass and of peak spring bloom amplitude were also estimated for each grid cell and each year. Using those indexes, we obtained the spatial pattern for the percentage of years with 1) higher spring bloom amplitude, 2) later onset of the spring bloom, 3) later onset of 245 stratification, 4) deeper mixed layer depth, 5) later occurrence of sufficient light conditions for building phytoplankton biomass, and 6) higher concentration of winter zooplankton biomass in response to tidal forcing (tidal scenario vs. non-tidal scenario).
Furthermore, we studied the changes in the spring bloom phenology in response to the spring-neap tidal cycle (i.e. whether spring or neap tide promote/hinder NPP). Considering that several spring-neap cycles may take place during the spring bloom 250 development, we studied the NPP difference during of spring bloom development between tidal scenario and M2 scenario in relation to the spring-neap tidal phase. The period of spring bloom development was defined as the time period with an increase in NPP from 12.5% to 87.5% of the maximum NPP. During this time period, we identified the occurrences (within a time window of one fortnight cycle) of positive/negative maxima of the NPP difference and temporally related the day of maximum difference to the adjacent day of the spring tide. This enabled us to evaluate the impact of spring-neap tidal cycles on spring 255 bloom phenology.

Quantification of limiting pattern of phytoplankton growth: light vs. nutrients
In ECOSMO, NPP is estimated as the sum of net primary production for all phytoplankton functional groups (Eq.2, denoted by j). For each functional group, the net primary production (NPP) is calculated by multiplying the maximum growth rate specified for the functional group ( ) with the minimum value ( ) of all limiting terms ( ) (Liebig's law, de Baar, 1994) and 260 the prevailing amount of phytoplankton biomass (standing stock, ) (Eq. 2). The limiting term ( ) for each growth resource is derived from the Monod equation (Monod, 1942), using the concentration of each growth resource( ) (Si: silicate-only for diatom growth, N: nitrogen, P: phosphorus, L: light) and the specific half-saturation constant (h) (Eq. 3). Further details of the nutrient limiting terms are given in Daewel and Schrum (2013). We hereafter call the minimum value of all limiting terms (Eq.4) the limiting value. The limiting value quantifies the availability of growth resources with a range of 0-1. The closer the 265 value is to 1, the more sufficient the resource is. Additionally, we identified the most limiting factors for each phytoplankton type ( ) (N, P, and L for flagellates; Si, N, P, and L for diatoms).
We analysed the limiting value to represent the environmental conditions of phytoplankton growth and the spatial and temporal dynamics of the most limiting factor.

Vertical distribution of phytoplankton: detection of subsurface maximum layer 275
The mixing intensity in the water column controls the distribution of phytoplankton and nutrients. As suggested by previous studies, phytoplankton may develop high subsurface concentrations in layers of low turbulence such as the pycnocline; production continues locally in low-turbulent zones as long as the growth requirements of nutrients and light are balanced (Cullen, 2015). In the stratified season, we differentiated the NPP generated in the surface layer (above 15 m) and in the subsurface layers, as a subsurface biomass maximum (SBM) emerged. The SBM was defined by its width, which was small 280 compared to the water depth, and was persistent in both time and space (Dekshenieks et al., 2001). In this study, we regarded layers deeper than 15 m as the subsurface. As an SBM necessarily includes local peaks, we first selected the depth at which the first-order derivative of biomass changed from positive to negative in the vertical biomass profile as a potential location for an SBM peak. To further identify the boundaries of the potential SBM, different strategies were applied depending on the number of vertical layers on either side of the potential SBM peak. If there were more than 5 vertical layers on either side of 285 the potential SBM peak, the vertical layer with the local maximum in the second-order derivative on each side of the potential SBM peak was recognized as the boundary of the SBM layer (Benoit-Bird et al., 2009). Otherwise, the adjacent layers were assumed to confine the potential SBM. The SBM peak could be no shallower than 20 m. We estimated the local background biomass value by linearly interpolating the biomass values of the upper and lower edges to the depth where the peak in biomass emerged. If the peak maximum biomass exceeded 1.5 times higher than the estimated background biomass in the respective 290 water column, the local vertical plankton biomass maximum was considered an SBM. Similar methods which has been applied to analyse phytoplankton (Chlorophyll a) vertical profiles in the German Bight and more details were laid out in Zhao et al.,(2019).

Identification of representative grid cells for spring-neap cycle impacts
In addition to tidal forcing, atmospheric forcing and bathymetry modulates stratification (e.g. Van Leeuwen et al., 2015) and 295 production pattern (Daewel and Schrum, 2017). Consequently, tidal impacts on stratification and primary production are subject to spatial-temporal variability. Furthermore, non-linear interactions among tidal constituents are pronounced in shallower waters, as suggested by Backhaus (1985) in inshore areas for the German Bight and Danish Coast. Although we preliminarily estimated the influence of the spring-neap tidal cycle via the difference in NPP between the tidal scenario and the M2 scenario, related responses would not necessarily be visible in a fortnightly cycle. To better associate the variation in 300 NPP with the spring-neap tidal cycle, we identified specific grid cells where both currents and biochemical factors displayed a distinguishable spring-neap cycle. Those locations were identified by using the estimated squared coherence between the power spectra (SCPS) of currents and NPP (Stoica et al., 2005;Welch, 1967). By adopting the SCPS method, we were able to select representative grid cells where both NPP and velocity showed obvious spring-neap cycles.

Spatial changes in mean production
The average annual NPP and the difference in NPP between the tidal and non-tidal scenarios are shown in Fig.2. The areaaveraged NPP increases slightly from 100.7 to 103.2 −2 −1 when tidal forcing is applied (Table 1); however, high spatial diversity in the sensitivity to tidal forcing is shown. Generally, estimated tidal impacts on NPP are highest in the stratified shallow North Sea, with a maximum response of up to 60 −2 −1 (Fig.2c). In the non-tidal scenario, high productivity is 310 restricted to the near-shore shallow regions along the British coast and the European continental coast (Fig.2a), which are the main regions where the euphotic zone reaches the bottom and nutrient remineralization fosters production throughout the year.
The primary production at the coast is additionally supported by estuarine-type baroclinic circulation in summer, which transports detritus and nutrient-rich bottom water towards the coast (Ebenhöh et al., 2004;Geyer and MacCready, 2014;Hofmeister et al., 2017). Tides cause a significant reduction in stratification in the shallow near-coastal areas of the North Sea 315 and in the EC at Dogger Bank and south of Dogger Bank and foster the development of tidal mixing fronts. Consequently, the production pattern changes notably when tidal forcing is considered. The primary production maximum is shifted further offshore towards the frontal region (Fig.2b). Large areas of the SNS, including Dogger Bank, the eastern BC, and the Danish Coast in the east, together with the NT, exhibit an increase in NPP when tidal forcing is prescribed. The shallow near-coastal areas in the south and the deeper areas in the NNS show a negative response of NPP to tidal forcing. A stronger negative 320 response is observed in the highly dynamic EC (Fig.2c). The NPP of Dogger Bank and the tidal mixing front area south and southeast of Dogger Bank responds the strongest to tidal forcing, with a mean change in NPP of up to 60 −2 −1 , nearly doubling local production. The amplitudes of the decreases in NPP in the negatively responding area are smaller than those of the increases in NPP, with amplitudes no more than 40 −2 −1 (Fig.2c); the largest amplitudes are in the EC. The intensity of this difference might be slightly sensitive to the consideration of inorganic SPM (see Appendix C). 325 The tidally induced change in NPP is associated with variations in the spatial distribution of the main limiting resources (limiting pattern) (Fig.3). Generally, in the tidal scenario, the area experiencing nutrient limitation decreases due to the enhanced mixing of inorganic nutrients into the euphotic zone, especially in the shallow North Sea where the bottom and surface mixed layer interact with each other. Simultaneously, light limitation increases. The predominantly light-limited regions, which are restricted to the shallow coastal regions in the non-tidal scenario (Fig.3a), expand offshore in the tidal 330 scenario (Fig.3b). Tidally induced resuspension and mixing of particulates and DOM into the euphotic zone result in dominant light limitation in almost the entire shallow North Sea (below 50 m depth) (Fig.3b). In contrast, in the surface layers of the stratified area, summer nutrient limitation is predominant, and the limiting value remains below 0.3 in both scenarios. The change from nutrient to light limitation in the SNS changes the limiting value to >0.4 in the tidal scenario, allowing better resource exploitation in these areas and sustaining NPP during summer. 335 The subdomain-division method described in section 2.4.1 identifies 7 different subdomains (Fig.4) that show characteristic responses to tidal forcing. Based on the division and the point-wise correlation of NPP variations in each subdomain (Fig.A1), representative grid cells were selected to study the mechanisms underlying the spatial variability of tidal responses in detail.
Areas with correlation coefficients higher than 0.3 occupied at least 53% of each subdomain, comprising 77% of the entire study area. This indicates that the division effectively explains the spatial diversity of the system with respect to the tidally 340 induced changes in NPP and the predominantly inherent similarity within each subdomain. The seven identified subdomains are listed below (Fig.4): 1. The English Channel (EC; dark blue) is characterized by an early onset of the spring bloom, strong mixing due to tidal stirring and shallow bathymetry. The EC is the most productive area in the non-tidal scenario (Fig.2a), with a mean NPP above 120 −2 −1 . 345 2. Negatively responding southern North Sea (neg. SNS; blue). The neg. SNS is separated from the EC by 52°N and from the positively responding area in the southern North Sea. The neg. SNS characterizes the permanently mixed area in the shallow water near the coast.
3. Positively responding southern North Sea (pos. SNS; light blue). This area includes the frontal regions that were identified as the areas with the highest responses in NPP (Fig.2). 350 4. Eastern British coast (BC; green). This area is a highly productive, positively responding inshore region of the eastern British coast. Sea (deep NNS; yellow). The deep NNS region coincides with areas of seasonal stratification and the lowest annual NPP in the tidal scenario (Fig.2). In this area, a slight decrease in NPP is estimated when tidal forcing is considered. 355 6. The Norwegian Trench (NT; orange) represents the area off the Norwegian coast, which is strongly impacted by the low saline outflow from the Baltic Sea. The NT shows a slight increase in NPP due to tidal forcing (Fig.2). 7. Low-sensitivity area in the northern North Sea (low-sen. NNS). The magnitude of the response of NPP to tidal forcing here is below 5 −2 −1 . This subdomain is influenced by two amphidromic points in the eastern North Sea, with tidal amplitudes of the M2 partial tide generally below 0.5 m. 360 Some narrow transient zones between the positively responding areas and negatively responding areas are shown in white in Fig.4. These transient zones with an absolute variation in NPP less than 5 −2 −1 are excluded from the following analyses.

Deeper northern North
Changes in NPP in response to tidal forcing for each subdomain are listed in Table 1.
The subdomain division corresponds well with the regional characteristics of M2 tidal energy dissipation rates, as suggested by the simulation study of Davies and Sauvel (1985). The EC subdomain includes the areas with the highest tidal energy 365 dissipation rates, which exceed 1000 −2 −1 (Davies and Sauvel, 1985). In most of the neg. SNS and some parts of the EC, the tidal energy dissipation rates are in the range of 100-1000 −2 −1 . In the pos. SNS, the BC and part of the deep NS, tidal energy dissipation rates range from 10 to 100 −2 −1 . The low-sen. NNS and NT are located in the area with tidal energy dissipation rates below 10 −2 −1 . The strong tidal energy in the SNS destabilizes stratification, as also revealed by the subdivision based on stratification patterns presented by Van Leeuwen et al. (2015). Our neg. SNS and EC subdomains coincide 370 with permanently mixed regions defined in the above study; in addition, the defined BC correlates with mixed or temporally stratified belts along the eastern British coast, as suggested by Van Leeuwen et al. (2015). The subdomains identified in the NNS coincide with seasonally stratified areas in the aforementioned study. However, the majority of pos. SNS, which shows the strongest response to tidal forcing, could not be identified with the method of Van Leeuwen et al. (Van Leeuwen et al., 2015). due to the variable stratification in these frontal areas induced by the spring-neap cycle, wind forcing, river runoff and 375 air temperature (Dippner, 1993;Schrum et al., 2003a;Sharples and Simpson, 1993). The subdomains also agree well with subdomains previously identified by ICES (ICES, 1983). Compared to the ICES subdivisions, which were determined considering biochemical and hydrographical characteristics (ICES, 1983), the northern 4 subdomains in our study coincide with regions where the gross water mass influx is mainly influenced by Atlantic water inflow. In contrast, for the 3 subdomains in the south, the influence of wind is more important for water mass exchange (Siegismund, 2001). 380

Characteristic seasonal changes
Out of the 7 subdomains (Fig.4), we selected 3 representative subdomains for further analysis of the changes in seasonality of NPP and the respective associated mechanisms. The neg. SNS represents the area along European continental coast where strong tidal forcing leads to permanent mixing and the NPP decreases as a consequence of tidal forcing. The pos. SNS embodies the transient zone between mixed and stratified water column and is characterized by the most significant positive response of 13 NPP to tidal forcing. The deep NNS is characterized by stable seasonal stratification. Here, the bottom mixed layer and surface mixed layer are well separated; thus, tides have a limited impact on the euphotic zone. The averaged time series  for each subdomain and the time series of the vertical profiles of each most representative grid cell (see section 2.4.1) are given in Fig.5 and Fig.6, respectively.
In the neg. SNS (Fig.5a), the spring bloom is delayed and strong fluctuations appear during the productive season in both 390 scenarios ( Fig.5a & Fig.6a, b, c, d). The pulses in NPP are probably due to predator-prey interactions and possibly modulated by advection. These pulses in NPP have previously been described by (Tett and Walne, 1995). The length of these fluctuations is slightly longer in the tidal scenario than in the non-tidal scenario, and changes in bloom initiation and the length of the quasi-periodic fluctuations generate positive-negative fluctuations in the NPP difference between both scenarios. We found no nutrient limitation in the water column in either scenario (Fig.6a, b) and no significant changes in the limiting values (Eq. 395 4) (note: the minimum limiting value stems from light limitation), except for the slightly higher values in deep water column under the tidal scenario (Fig.6. c, d). This exception is likely caused by the downward mixing of shade-producing organic materials (e.g., phytoplankton, DOM and detritus), which leads to improved light conditions in the upper layer and better penetration. However, this result does not explain the negative NPP response in the area. Lower NPP in the tidal scenario than in the non-tidal scenario, especially in spring and early summer, results in an overall negative response in NPP. A likely reason 400 for the reduction in NPP in the neg. SNS subdomain could be the tidally induced dilution of phytoplankton biomass in the euphotic zone in the shallow areas. The increased mixing in the tidal scenario dilutes the phytoplankton concentration in the upper, highly productive water layer (see vertical profiles of biomass Fig A4, a,b) and consequently reduces the time during which phytoplankton cells are exposed to high surface irradiance. Considering the small difference in the growth resources between the two scenarios (Fig.6c, d), we mainly attribute the variation in NPP to the vertical distribution of standing stocks. 405 The most dominant change in seasonality as a consequence of tidal forcing in the seasonally stratified subdomains (pos. SNS and deep NNS) is the delay of the spring bloom in the tidal scenario (Fig.5b, c). However, in the pos. SNS, this delay is only a few days long; in the deep NNS, this delay encompasses one month. Accompanying the delay, the amplitude of the spring bloom in the tidal scenario, especially in the pos. SNS, exceeds that of the non-tidal scenario. The spring bloom in the NS typically consists of diatoms, while after silicate depletion, flagellates dominate the summer production (McQuatters-Gollop 410 et al., 2007;Schrum et al., 2006). Comparing the seasonality of NPP variation (Fig.5b, c) with the annual averaged NPP deviation between scenarios (Fig.2c), we found that the variation in NPP in summer is basically in phase with the direction of the NPP's response to tidal forcing for both the deep NNS and pos. SNS. Especially in the pos. SNS (Fig.5b), summer blooms are higher in the tidal scenario than in the non-tidal scenario, with a maximum difference in July and August, fostered by weaker stratification and regular nutrient injections into the surface mixed layer due to tidally induced turbulence (Fig.6f, h). 415 Surface summer production is sustained throughout the summer at values of approximately 50 −3 −1 and more in the upper 15 m (Fig.6f), and light remains the dominant limiting factor in the surface layer, except for a temporal silicate limitation after the spring bloom (Fig.6h). In contrast, without tidal stirring, surface waters become nutrient depleted soon after the spring bloom in May. After silicate limitation, nitrogen limitation persists (Fig.6g) in the surface waters throughout the seasonal and enhanced turbidity caused by tides, no SBM production occurs in this simulation. The nutrient supply advantage in the tidal scenario persists until the beginning of October (Fig.6f), when the water column in the non-tidal scenario is also mixed by atmospheric conditions, causing an increase in production at the surface. For the pos. SNS the modulation of nutrient availability is the most important factor responsible for changes in NPP. The high biomass stays in the pycnocline during summer due to the weak mixing in the non-tidal scenario. In contrast, due to the weak stratification and strong mixing, 425 generated high biomass is continuously mixed in the euphotic zone in the tidal scenario (Fig.A4, c,d).
In the deep NNS, the influence of tides on NPP is relatively weak and mainly visible in summer (Fig.5c). The deep NNS (Fig.6i-l) is typically characterized by stable seasonal stratification and summer subsurface primary production in both the tidal and non-tidal scenarios. The delay of the spring bloom in the tidal scenario causes a quicker succession and consequently overlapping diatom and flagellate blooms (Fig.6l, Fig.5c). The productive period, which lasts nearly 3 months and includes 430 two pulses of NPP in the non-tidal scenario (Fig.6i), is shortened to 6 weeks in the tidal scenario (Fig.6j). The NPP contributed from subsurface production is higher in the tidal scenario than that in the non-tidal scenario (Fig.6i, j). Because of stratification and nutrient depletion, high biomass is confined to a region within the pycnocline in both scenarios. The SBM in the tidal scenario deepens because of mixed layer deepening due to tides (Fig.A4, e,f).
In the other identified subdomains (results not shown), the changes in primary production basically follow the pattern explained 435 above. In the EC subdomain, the tidal impact on production is comparable to that in the neg. SNS, whereas in the BC subdomain, nutrients are rarely the most limiting factors due to weak stratification, and the response can be compared to that in the pos. SNS. In the low-sen. NNS, where tidal dissipation is weak, the vertical distribution pattern of NPP in both scenarios is almost identical.
Our results indicate that, in principle, tidal stirring causes two major changes in the NPP pattern: i) a change in the spring 440 bloom phenology of some areas and ii) an altered ratio between surface and subsurface production. Both features merit further discussion, which is given in the following paragraphs.

Changes in spring bloom phenology
As one of the most important biological events in the NPP annual cycle (Bagniewski et al., 2011;Sabine et al., 2004), the spring bloom requires specific attention. As shown by the time series analysis for some subdomains (Fig.5) and the time series 445 of profiles at the representative points (Fig.6), the postponement of the spring bloom is a prevalent phenomenon when tidal forcing is applied. The changes in spring bloom phenology and the processes responsible for these changes , such as the delay in the onset of stratification, variations in light conditions, the mixed layer depth and winter zooplankton concentrations (Fig.7), were analysed using the method outlined in section 2.4.2.
In line with the distribution of tidal energy dissipation given by Davies and Sauvel (1985), the spring bloom delay is robust in 450 the SNS and along the British coast (Fig.7a), while in the northeastern part of the North Sea, the spring bloom is delayed in no more than 50% of all years. An increase in the peak spring bloom biomass (Fig.7b), mainly in areas with a positive response of NPP to tidal forcing (Fig.2c). However, in some isolated locations in the negatively responding areas, such as the neg. SNS and EC, the spring bloom amplitudes are still higher in the tidal scenario than that in the non-tidal scenario in more than 50% of the years. One potential reason for the spring bloom delay is a change in light conditions, especially in very shallow coastal, 455 non-stratified areas where tidal stirring enhances resuspension in the water column (Fig.7c). The onset of light conditions sufficient for phytoplankton growth in the well-mixed water column is delayed in the coastal areas of the southern and eastern boundary and in the shallower parts of Dogger Bank. However, the distribution of this impact does not explain the major patterns of changes in the spring bloom phenology. Tides also increase mixing and hence potentially prevent stratification in shallow water columns or delay the onset of stratification, as discussed previously by a number of authors (Bowden and 460 Hamilton, 1975;Loder and Greenberg, 1986). Because tidally induced energy dissipation is cubically proportional to the strength of tidal currents (Simpson and Hunter, 1974), we can expect the strongest variation in stratification in regions with the strongest tidal currents, as observed along the British coast, in the EC and in the German Bight (Davies and Sauvel, 1985). This expectation is supported by earlier observations suggesting that the onset of the spring bloom is triggered by improved light conditions because of solar radiation and stratification (van der Woerd et al., 2011). The onset of stratification (Fig.7e) 465 in the tidal scenario is mainly delayed in the Scottish coast and the frontal areas of the SNS. Furthermore, the response of stratification to tidal forcing is more stable in the southwestern part (the Estuary of Humber, Dogger Bank) than in the southeastern part of the SNS (Fig.7e). Apart from solar heating, the stratification in the southeastern part of SNSis additionally influenced by freshwater supplies from land and wind forcing (Jacobs, 2004;Ruddick et al., 1995;Schrum, 1997).
Consequently, the variation in the onset of stratification is less clear in the southeastern part than in other parts of the SNS. In 470 the NNS, the tidal wave propagation deepens the mixed layer depth (Fig.7d), which similarly results in a later onset of the spring bloom, despite only weak changes in the onset of stratification. As a consequence of the thicker layer in which phytoplankton are mixed, the phytoplankton are less exposed to the favourable surface light conditions and will thus take longer to build up the spring bloom biomass.
Although the North Sea is in principle a bottom-up-controlled ecosystem, zooplankton predation is occasionally an important 475 process controlling NPP (Daewel et al., 2014). In early spring, even under favourable growth conditions, the spring bloom will only initiate until production exceeds the loss due to grazing (George et al., 2015;Martin, 1965). This grazing pressure is basically correlated with the overwintering zooplankton stock. Based on our results, increases in the winter zooplankton biomass and delays in the spring bloom coincide only in the frontal region of the SNS and central NNS. Therefore, we conclude that the delay in spring bloom by tides is mostly due to bottom-up control. 480 The spatial pattern given in Fig.7 shows that the delayed onset of the spring bloom in the tidal scenario may mainly be attributed to deteriorated light conditions in the shallow well-mixed area (Fig.7c) and changes in the stratification of seasonally stratified areas, such as delays in the development of stratification (Fig.7e) or the deepening of the upper mixed layer (Fig.7d).
Although the predator biomasses are higher prior to spring bloom in some areas, enhanced grazing pressure at the beginning of the bloom period does not seem to be the main mechanism delaying the onset of the spring bloom (Fig.7f), although we 485 assume this pressure plays an additional role in the central NNS and frontal regions.

Changes in subsurface production in stratified season
To further quantify the magnitude of the changes in surface and subsurface production during the stratified season, we separated NPP vertically into upper-layer production (above 15 m) and production in the SBM layer and compared the results between scenarios (Fig.8), using the mean annual value for the analysed period . At the stratified side of the frontal 490 zones (pos. SNS), the surface production response of NPP is positive almost everywhere, with a maximum reaching +50 −2 −1 (Fig.8b) at south of the Dogger Bank. In contrast, the changes in response to tidal forcing within the SBM show both negative and positive responses around the Dogger Bank (Fig.8a). A positive response to tidal forcing, which is generally one order of magnitude smaller than the increased amplitude of NPP in the surface layer, occurs only at the northern edge around Dogger Bank and the deeper part of the German Bight. A similar pattern with a strong positive response to tidal forcing 495 at the surface and a negative response in the SBM appears in the BC area. In line with former studies in the North Sea, the NPP in the upper layer dominates the whole production budget (van Leeuwen et al., 2013). Although the expansion and duration times of the SBM decrease due to tidal forcing, e.g., in the inshore areas along the BC and at the Danish Coast ( Fig.A3,b,d), tidal forcing promotes NPP within the SBM in some areas, especially at the northern edge of the Dogger Bank.
Observational studies suggested that the productive areas at the edge of the Dogger Bank are fueled by baroclinic circulation 500 related to the front and the spring-neap adjustment (Pedersen, 1994). When considering an SBM duration of 110 days (Fig.A3, c) at the northern edge of the Dogger Bank, the average daily NPP (deduced from the annual NPP (Fig.A3,a)) is approximately 239 −2 −1 , which corroborates the observation-based estimate of NPP (295 −2 −1 ) calculated from measured oxygen surplus concentration data (Richardson et al., 2000).
In the NNS, the variation caused by tidal forcing in NPP is below 15 −2 −1 (Fig.2c). In some parts of the deep NNS, the 505 tidal forcing causes higher production in the SBM and lower production at the surface (Fig.6i, j; Fig.8a). Due to the decoupling between the surface and bottom mixed layers, the pycnocline acts as a barrier that keeps the stirred-up nutrients below the pycnocline and sustains NPP in the SBM (Fig.6i, k). Because the amplitude of NPP variations in the upper layers is ten times higher than that in the SBM (Fig.8), the overall response to tidal forcing is negative (Fig.2c) in the deep NNS.

Impacts of the spring-neap cycle 510
The spring-neap tidal cycle introduces a fortnightly periodic change in tidal mixing, which has a significant influence along the British coast and in the English Channel (Fig.9). The differences in current speed between the tidal and M2 tidal scenarios vary over the spring-neap tidal cycle. The maximum spring-neap range of these differences is up to 0.3-0.6 m/s (Fig.9), indicating that a non-negligible change in turbulent kinetic energy is introduced to the water column via the spring-neap cycle.
Here, we will provide model estimates on the spatial variability in the resulting response of the NPP to the spring-neap cycle 515 and explore the potential mechanisms of these responses.
Annual NPP changes induced by the spring-neap cycle reach maximum values of up to 5 −2 −1 (Fig.9). Although this amount is relatively small compared to the overall system productivity, the changes due to spring-neap dynamics could be very relevant locally and in specific time periods. An average positive response of NPP emerges in the southeastern part of the North Sea, in the English Channel and along the British Coast (Fig.9). The highest mean changes in NPP are found in the 520 western part of the Dogger Bank, in the English Channel and off the Scottish coast. In contrast, a negative response in annual production emerges off the Northumbrian coast and in the Southern Bight off the European continent (Fig.9). The response of NPP to spring-neap tidal forcing is weak in early spring and winter (data not shown). Under mixed conditions or during periods of the establishment and decay of stratification, spring-neap tidal mixing can be overridden periodically by other mixing events (e.g. driven by wind); hence, pronounced irregularities in NPP responses to spring-neap tidal forcing are detected. A significant 525 response of NPP to spring-neap tidal forcing is found for summer periods under stable stratification. To illustrate the basic mechanisms responsible for the response of NPP due to spring-neap tidal cycle, we present time series of the biomass, nitrate, NPP and turbidity (Eq. 1) profiles for two characteristic grid cells (selection described c.f. 2.4.5) that respond differently to spring-neap tidal forcing. The near-shore grid cell off the Estuary of Humber (EH, Fig.9) shows a negative response, and a grid cell located at the frontal zone at the western edge of Dogger Bank (WDB, Fig.9) responds positively to spring-neap tidal 530 forcing. The model results are presented for a couple of selected successive spring-neap tidal cycles simulated for the year 2001 (Fig.10).
The EH site, which is located further inshore compared to the WDB, is characterized by high turbidity. The increased nitrogen in the upper layers is in phase with elevated turbidity but in anti-phase with biomass and NPP. This phenomenon indicates that during spring tide, the process of phytoplankton biomass dilution (Fig.10c) and shading due to the upward mixing of organic 535 material (Fig.10b) slows NPP in the upper mixed layer, resulting in a negative NPP response during spring phases (Fig.10d).
The elevated NPP reaches a maximum at the end of the neap phase (Fig. 10d), possibly because of the reduced vertical mixing.
The decreasing turbidity in the neap phases, despite increases in phytoplankton biomass, reveals that suspended and resuspended organic material have a reduced impact on the surface light conditions during neap phases compared to the spring phases (Fig.10b). In neap tidal phase, given less vertical mixing, phytoplankton cells remain in the lighted surface layer for 540 longer time and access better light conditions; hence, the available nutrients can be utilized for phytoplankton growth (Fig.10c).
In contrast, the WDB site is typically characterized by seasonal stratification and summer nutrient (i.e., nitrate) depletion in the surface layer. However, as the WBD site is located in the frontal zone, relevant factors in this zone do not necessarily show the spring-neap fluctuation as clearly as those in the EH site. During spring tide, enhanced vertical mixing dilutes the phytoplankton biomass in the upper layer and redistribute biomass more evenly in the whole water column, resulting in less 545 phytoplankton biomass in the upper layer (blue) and more biomass in the lower layer (red) compared to that in the M2 scenario (Figgure.10g). Spring tidal forcing results in the replenishment of nutrients in the euphotic zone and a pulse of increased NPP follows spring tide mixing (Fig.10e, h). The downward mixing of biomass into lower layers has no substantial negative effect on NPP during spring tide (Fig.10g). As a consequence of nutrient replenishment in the surface layer during former spring tides, given less vertical mixing during neap tide, biomass increases in the upper layers (Fig.10g). Resuspension effects 550 resulting in increased turbidity at lower layers are visible from neap to spring but do not significantly change turbidity in the surface layers (Fig.10f). Surface turbidity changes are consequences of increased NPP (Fig.10f).
Observation based estimates of spring-neap impacts on NPP given by Richardson et al. (2000) found that increased nitrate fluxes by tidal pumping contributed to NPP with 4-6 −2 for one spring-neap cycle at the northern edge of the Dogger Bank, mainly due to increased production in the subsurface layer. By upscaling these results to the entire stratified season, 555 considering 6-8 spring neap cycles, Richardson et al. (2000) proposed that the additional NPP contribution by the spring neap cycle was in the range of 24-48 −2 for the whole stratified season. We resampled the simulated NPP along the same transect as sampled by Richardson et al. (2000) and for the same time period (29/07-04/08 1997). We extended the time period to 26/07/1997-08/08/1997 to cover a full spring-neap cycle and found our simulated response of NPP to tidal forcing (the tidal scenarionon-tidal scenario) is 3.03 −2 for one spring-neap cycle (Fig.A5,a). These values are slightly below the lower 560 edge of Richardson's estimates (4-6 −2 ). However, simulated frontal locations are not always conform to the observed fronts due to unresolved sub-scale processes, which remain unconsidered in a 10km x 10 km model resolution and coarse atmospheric forcing (NCEP/NCAR reanalysis). When we resampled the NPP along the fronts in our simulation, which is at a distance of a few grid points further south from the fronts in Richardson et al., (2000) (Fig.A5,a), we found that the simulated change in NPP (5.99 −2 for one spring-neap cycle) reaches the upper level of estimates based on observations (Fig.A5,a, 565 Table.A1). When we compare the NPP response throughout the whole stratified season (simulated as 15 −2 ) we find this to be lower than Richardson's upscaling estimation (24-48 −2 for the whole stratified season). The reason for this discrepancy is a too simplified upscaling procedure used by Richardson et al., (2000), neglecting the sensitivity to seasonality.
Conditions measured over a few days between July and August (Richardson et al., 2000) are not representative for the whole stratified season. In contrast to Richardson et al.'s conclusion (2000) that spring-neap cycle played the major role in fuelling 570 NPP, our study indicates further, that the semidiurnal tide plays the major role in pumping up nutrients and sustaining the NPP but not the spring-neap cycle as hypothesized by Richardson et al. (2000). Our estimates of averagely 0.14 −2 of NPP promoted by the spring-neap tidal cycle during one tidal cycle (Fig.A5,b) is considerably lower than that supported by the standard tidal forcing (M2+S2) (5.99 −2 ) (Fig.A5,a) and hence spring-neap tidal pumping contributes only little to the increase in NPP. Based on our simulation, tidal pumping sustaining subsurface NPP mainly occurs in July and August, with 575 an average value of approximately 3 −2 ℎ −1 in frontal areas around Dogger Bank. This is close to the estimate of Richardsons et al. (2000). In other weakly stratified months, the value is no more than 1 −2 ℎ −1 or even negative (data not shown). Sharples (2008) investigated a similar question for the Celtic Sea with model simulations. He found that the NPP varied up to 580 70% with the spring-neap tidal cycle. We could not confirm such high tidal impacts on NPP for the North Sea; our estimates of the response of NPP to the spring-neap tidal cycle are only up to approximately 10% of the tidal impact (M2+S2) on NPP ( Fig.2 & Fig.9). One explanation for this discrepancy is the higher spring-neap tidal current amplitude in the Celtic Sea compared to the North Sea, which may result in a stronger response of NPP to the spring-neap cycle. However, it is also possible that the simpler model setup used by Sharples, such as neglection of advection, a constant grazing rate and neglected 585 impacts on resuspension and shading by DOM and detritus, resulted in higher NPP sensitivity to tidal forcing in their simulation.
As discussed in section 3.2, tidal forcing not only impacts the magnitude of NPP but also spring bloom phenology. It is reasonable to assume that spring-neap tidal forcing also modulates the development of the spring bloom. To understand the impact of the spring-neap phase on the biomass build-up during the spring bloom, which typically occurs over one or several 590 spring-neap cycles, we related the time periods with the maximum difference in NPP between the tidal scenario and the M2 scenario to the spring-neap cycle phase (Fig.11) at the SN site (see Fig.9). The SN site is located in the tidally energetic northwestern North Sea, where the development of the spring bloom often benefits from thermal stratification (Rodhe, 1998) but is sensitive to episodic 'noise' added by wind forcing (Waniek, 2003). During the development of the spring bloom, in the difference between NPP time series (tidal scenario-M2 scenario), an increase in NPP often occurs in neap phases, whereas 595 NPP is often decreased in spring phases. This indicates that the development of the spring bloom benefits from the neap phase but is interrupted or dampened during the spring tide (Fig.11). A similar phenomenon has been explored and confirmed by Sharples et al., (2006) at a site south of the SN site. As suggested by Sharples et al., (2006), the onset time of the spring bloom is shifted by the spring-neap tidal cycle because the onset or intensity of stratification is strengthened during neap tides when the vertical mixing is dampened. 600

Summary and Conclusions
A model-based sensitivity experiment with varied tidal forcing was performed to evaluate tidal impacts on NPP, considering the major bottom-up controlling processes, including the tidal mixing of nutrients, organic matter and plankton biomass and tidal resuspension of suspended matter. The responses to tides in the North Sea differ regionally and depend on the local hydrodynamic characteristics. In permanently mixed areas in the southern part of the North Sea, light availability is the major 605 limiting factor. The enhanced tidal resuspension and mixing of suspended matter into the surface layers deteriorate light conditions in the upper layers for phytoplankton growth and thus hinder primary production. In contrast, in frontal areas and seasonally stratified areas in the SNS where stratification is susceptible to tidal mixing, nutrient replenishment due to tidal forcing sustains NPP in summer and thus contributes a significant increase in NPP in both the surface layer and within the pycnocline. In the NNS, which is characterized by relatively weak tidal forcing and deep bathymetry, the bottom and upper 610 mixed layers are well separated, and the influence of tidal forcing on NPP is limited.
However, the quantitative estimates provided here are model-and parameterization-specific. Dominant biochemical processes are generally well represented in simplified NPZD-type models, and the ECOSMO model used here is applicable for resolving ecosystem dynamics at seasonal to decadal time scales when forced by realistic boundary conditions (Daewel and Schrum, 2017). However, parameterization and unconsidered processes, such as the role of macrobenthos in the system, internal waves 615 at the shelf break, and coastal light attenuation due to inorganic suspended matter, and simplified physiological processes could potentially modulate or change the model's sensitivity to tidal forcing. Studies identifying the contribution of these processes to tidal impacts on primary production are needed; thus far, we can only speculate on potential impacts.
Macrobenthic grazing likely changes the biochemical cycling and turbidity in the water column, subsequently changing the sensitivity of NPP to tidal forcing. In shallower waters, high near-bottom concentrations of suspended organic matter are 620 susceptible to mixing into the euphotic zone, and increasing light attenuation leads to decreasing production (see Fig.7c, cf. section 3.1). Macrobenthic biomass, specifically from filter feeders, might significantly reduce resuspension and near-bottom suspended matter concentrations, thereby increasing the proportion of organic matter that remains in the food web (Prins et al., 1996). From observations, we know that macrobenthos show a distinct spatial pattern following principle production patterns in the North Sea with higher biomass in the shallow SNS (Heip et al., 1992). Therefore, we can expect an increase in 625 NPP sensitivity to tidal forcing due to macrobenthos activity.
Furthermore, we hypothesize that the positive response of NPP to tidal forcing in the NNS was underestimated by our simulation due to the implementation of identical boundary conditions for all scenarios. We neglected the influence of tidalgenerated internal waves on nutrient conditions. Tidal-generated internal waves are initiated at the shelf edge and enhance turbulent mixing at the shelf break and on the shelf (Heathershaw et al., 1987;Loder et al., 1992;New and Da Silva, 2002;630 Sharples et al., 2001). As internal tides break at the shelf edge, energy dissipates mainly at the shelf break and other bathymetric features, which causes vertical mixing that drives vertical nutrient fluxes and sustains phytoplankton growth (Holligan et al., 1985;Pingree et al., 1981;Sharples et al., 2007). Therefore, internal tidal waves will likely lead to mixing and increase nutrient pulses onto the shelf, consequently supporting NPP. In our setup, the average impact of tidal-generated internal waves on nutrient concentrations was considered with the climatological boundary conditions (Conkright et al., 2002), and differences 635 among the simulated scenarios were not considered.
Another source of uncertainty in our model stems from the neglection feedbacks related to inorganic material, which influences underwater light conditions, especially in shallow areas. Seasonal differences in yellow substance concentrations coincide with fresh water input (Schaub and Gieskes, 1991;Warnock et al., 1999). There are two main sources of SPM plumes in the North Sea. One source lies at the southern British coast and originates from local discharges (Humber-Wash & Thames rivers), 640 coastal erosion and influx from the English Channel (Eisma, 2009). The other major source of SPM originates from the large continental rivers and diffusive sources entering the North Sea from the European continental coast, particularly off the Belgian coast and the Wadden Sea (van Alphen, 1990;Postma, 1981). Waves and currents are the controlling factors of the dispersion, resuspension and deposition processes of SPM (Holt and James, 1999). In winter, the two SPM plumes expand further offshore due to intensified mixing and both SPM plume deposits in both the Skagerrak and Norwegian Channels. We have evaluated 645 the potential impacts of inorganic SPM on our findings (Appendix C) and found that the general results regarding the tidal impacts on NPP remains largely insensitive to consideration of inorganic SPM and its seasonality. The reason lies in the spatial and temporal distribution of the SPM in the North Sea. During summer, SPM concentrations are low (Fig.C1)especially in stratified conditions and upper water layers (Capuzzo et al., 2013;Dobrynin et al., 2010). Only in shallow areas with permanent mixing concentrations are higher (Van Raaphorst et al., 1998). This is critical for our analysis since most differences in NPP actually occur in summer stratified conditions. A simulation study (Tian et al., 2009) in the German Bight found that implementing SPM is only critical at the onset of bloom, given reasonable parameterization, similar bloom amplitude was achieved in scenarios including or omitting SPM. Furthermore, measurements suggested that in the central North Sea, the water body itself triggers most of the attenuation (Jones et al., 1998). SPM is more relevant to attenuation in near shore areas due to cliff erosion and river input (Eisma, 2009). The relevance to turbidity of fluvial SPM is confined to river mouths because 655 SPM deposits quickly (Pleskachevsky et al., 2011;Siegel et al., 2009). Organic suspended matter (which is considered in the model) accounts for a high fraction of the total suspended matter (TSM) in most areas in the southern North Sea except for the very near shore areas (Schartau et al., 2018). The areas were inorganic suspended matter dominates are in the negative responding regions of our analysis (Fig.2c). The distribution of inorganic suspended matter is influenced by many factors, such as transportation with residual currents, aggregation with organic matter, type of benthic sediments and so on. Clearly, 660 interaction processes as mentioned above cannot be resolved by implementing a climatological SPM field. Thus, the numerical experiment presented here is considered to be a first step towards understanding the role of SPM for tidal impacts, and further studies specifically focussing on shallow coastal areas would require reasonable boundary conditions for inorganic matter from benthic sediments and river inputs as well as a more reasonable representation of bio-physical interactions related to inorganic matter. However, this is beyond the scope of the current study and should be emphasized more thoroughly in future work. 665 The major tidal impacts on NPP are via vertical mixing. Given the small horizontal gradient of both nutrients and biomass and weak tidal residuals of no more than a few centimetres per second (Prandle, 1984), the impacts of horizontal advection are negligible. To investigate the influence of advection on the concentration of nutrients and phytoplankton biomass in our study, we estimated the net horizontal transport between grid cells. In most parts of the study domain, we found that the contribution of mean tidal advection does not exceed 5% (not shown here). Exceptions occur in the Skagerrak Channel, where relatively 670 high residual currents drive water exchange between the North Sea and Baltic Sea (Brettschneider, 1967), and in the EC close to the model boundary, where relatively high current speeds caused by atmospheric forcing and topography emerge irregularly, mainly in spring and winter. However, this result is not true for smaller horizontal and temporal resolutions.
Since the North Sea can in general be considered as bottom-up controlled (Daewel et al., 2014;Heath, 2005), using a lower trophic level model for investigating tidal impacts on NPP is a valid approach. Although situations with clear top-down control 675 on zooplankton have been observed (Munk and Nielsen, 1994), these events occurred highly restricted in time and space and assumed to be only of minor relevance for the general processes described in this manuscript. In previous studies, which addressed similar scientific question, constant grazing rates (Sharples, 2008) or grazing loss being proportional to phytoplankton biomass (Cloern, 1991) were prescribed in the simulations. In this study, we utilize a lower trophic level NPZDtype model, only considering lower trophic level dynamics up to zooplankton, which is simulated as a state variable considering 680 feeding preference, growth, excretion and mortality. Fish predation is only implicitly considered as part of the zooplankton mortality rate. Simulations with ECOSMO E2E (an updated version of the ECOSMO model) including functional groups for fish and macrobenthos revealed that temporal and spatial variations in zooplankton mortality due to fish predation are determined by the specific hydrodynamics of the North Sea (Daewel et al., 2018). Repeating a similar study with an NPZD-Fish model would be interesting, however, beyond the scope of our study. 685 Given the importance of tidal forcing to NPP, especially in frontal areas, which are known to be biological hotspots (Belkin et al., 2009), tidal impacts on higher trophic levels than those studied here merit further consideration and investigation in the future. Regarding the growth of macrobenthos, tidal stirring influences the sinking and resuspension of organic matter and thus influences food quality and bioturbation (Foshtomi et al., 2015;Zhang and Wirtz, 2017). Tidal forcing in frontal areas not only provides enough prey for fish larvae due to nutrient enrichment and higher NPP but also influences convergence zones, 690 which are typical places for fish spawning and nursing (Bakun, 2006). Further investigations based on a combination of observations and multi-process coupled simulations could enable a better understanding of the impacts of tidal forcing on ecosystem processes and their variability. Long-term tidal variations, such as the 18.61-year nodal cycle or the 8.85-year lunar perigee cycle, merit particular consideration.      Contour lines indicate the estimated mean spring-neap cycle range of the tidal current speed difference between the tidal (M2 + S2) and M2 scenarios. The two magenta dots indicate the locations of 2 characteristic grid cells. One grid cell is close to the Estuary of Humber (EH), and the other grid cell is located more offshore at the western edge of Dogger Bank (WDB) (see Fig. 10). The magenta star shows the location of the grid cell used for the analysis of the advancement and delay of the spring bloom due to spring-neap tidal forcing (see Fig. 11). Figure 10. Spring-neap cycle impact on nitrate (a, e), turbidity (b, f), phytoplankton biomass (c, g) and primary production (d, h). Differences between the two scenarios (tidal scenario (M2+S2) -M2 tidal scenario) are presented for two characteristic points, i.e., EH: upper panel (a, b, c, d), WDB: lower panel (e, f, g, h). To show the periodical fluctuation of currents and NPP and relate these fluctuations to changes in nitrate, turbidity, biomass and NPP, the differences in the depth-averaged velocity amplitude (black) and depth-integrated NPP (green) are presented in each subplot; both time series underwent smoothing with a 24-hour running mean. Figure 11. The occurrence of an increase (red) or decrease (black) in the NPP difference (tidal scenario (M2+S2)-M2 tidal scenario) relative 1110 to the nearest spring tide (spring and neap phase indicated). The development of the spring bloom period is defined as the time when NPP increases from 12.5% to 87.5% of the maximum NPP prior to the major peak of the spring bloom.    in the simulation, we implemented a climatological SPM field for the North Sea (available with daily resolution and 31 vertical layers) for our simulation. This SPM field was derived from statistical regression model which considers tidal currents, salinity and water depth (Heath et al., 2002). The SPM field is able to resolve spatial distribution pattern and seasonal cycling of SPM concentration in the North Sea ( Figure. C1) and has been applied in many hydrodynamic-1280 biogeochemical coupled models (Große et al., 2016;Kerimoglu et al., 2017).
Taking the parameterization scheme proposed by Tian et al. (2009), we parameterize shading effects due to SPM as: = • √ (Eq.C1) 1285 The was set as 0.02 2 −1 . We added the contribution of SPM to the light shading scheme as described in the paper (Eq.1). We decreased the background attenuation co-efficient 1 , 0.03 −1 , to 0.025 −1 ( 2 ) and use the following parameterization for light attenuation: We implemented the new light shading scheme (Eq.C2) and evaluated the difference in NPP contributed by tide, by comparing the annual mean NPP in tidal and non-tidal scenarios using the Eq. C2 ( Figure C2). The general pattern remains largely insensitive to consideration of spatial and seasonal variations in SPM. The positive and negative responding areas hold the same distribution pattern, but the NPP's increasing amplitude with tidal forcing in frontal areas decrease slightly when 1295 SPM is explicitly considered. This is because the elevated NPP fueled by pumped up nutrients is partly offset by increased shading effects due to SPM. However, the sensitivity to SPM is minor and do not affect the general results of our study. both with SPM filed implemented. The difference in the mean annual NPP of both scenarios is in (c). The spatial coverage is smaller than original simulation domain since the SPM filed data is available from 50.5-57.5°.