Impacts of the Assimilation of Satellite Sea Surface Temperature Data on Volume and Heat Budget Estimates for the North Sea

Mechanisms controlling the heat budget of the North Sea are investigated based on a combination of satellite sea surface temperature measurements and numerical model simulations. Lateral heat fluxes across the shelf edge and into the Baltic Sea as well as vertical ocean-atmosphere heat exchange are considered. A 3-D variational (3DVAR) data assimilation (DA) scheme is applied, which contains assumed model error correlations that depend on the mixed layer depth derived from a coupled circulation/ocean wave model. The analysis balances pressure gradients introduced by temperature modifications. Significant hydrodynamic model response to DA was found, which should be considered in the heat budget estimations. The observed change of the current velocity field decreases the lateral advective volume/heat exchanges between the North Sea and the Atlantic, yielding an increased heat flux from the Atlantic into the North Sea and more heat flux from the sea to the atmosphere. The largest DA impact on volume/heat transport is in the Norwegian Channel, where the dominant process is Eulerian transport, followed by tidal pumping and wind pumping. Further analysis reveals an acceleration of the along-shelf current at the northern edge of the North Sea, a decrease in the horizontal pressure gradient from the Atlantic to the North Sea, and a reduction of the Eulerian transport of volume/heat outward the North Sea. Furthermore, the coupling between the circulation model and the wave model has significant impacts on lateral heat advection in the DA run, which is due to the wave impact on the mixed layer depth. the a data to of on and the of a and to investigate the impacts of sea surface (SST) on the and the We found that DA improved SST modeling, thereby modifying the volume and heat budgets between the North and the Atlantic. The largest change occurs at the Norwegian Channel, where the total water/heat transport from the North Sea outward is reduced. Moreover, SST assimilation also changes the air-sea heat This study improves our understanding of the relations between and

The effects of anthropogenic warming are known to be superimposed with natural variability patterns, such as the Atlantic Multi-decadal Oscillation (AMO) (Knight et al., 2005). To better understand the temperature dynamics in the North Sea, it is important to consider different components of the coupled system, such as heat fluxes between the atmosphere and the ocean, and the lateral advection of heat across the North West Shelf (Schrum & Backhaus, 1999). However, the exchanges that occur across the shelf edge are complex and the subject of ongoing research (Huthnance et al., 2009). Moreover, because of the strong spatial variations in warming patterns, the estimation techniques for regional SSTs must be further improved. One valuable source of information is passive microwave observations obtained from satellites, which provide SST maps on a kilometer scale. However, limitations are observed regarding the temporal and spatial effects of clouds coverage.
Assimilating data into numerical models and combining available information about the underlying physics with observations into an optimal estimate represents an efficient approach. In the pioneering study by Annan and Hargreaves (1999), a simplified Kalman filtering approach was applied to a three-dimensional (3-D) baroclinic model of the North Sea using satellite SST data. By using an anisotropic recursive filter scheme, Liu et al. (2009) assimilated the temperature and salinity profiles in a coastal ocean model and achieved a considerable improvement in the oceanic forecasting efficiency and accuracy. Fu et al. (2011) assimilated the temperature and salinity profiles in a regional ocean model but implemented an ensemble optimal interpolation. Other studies investigated the impacts of the timing and frequency of data assimilation (DA) (Losa et al., 2012) and the uncertainties of initial states (Losa et al., 2014) on the SST prediction performance. Many of these early studies focused on the predictive skills of models or the mathematical/statistical aspects of the assimilation technique. Furthermore, standard metrics such as the root mean square error (RMSE) have been used for prognostic model variables with respect to observations to assess the analysis performance. In contrast, the aim of the present study is to analyze the impacts of DA specifically on the water volume and heat budget in the North Sea. The following questions are at the center of our investigation: • How does the assimilation of SST observations change the different components of the simulated North Sea heat budget? • What are the secondary effects of temperature analysis on the remaining prognostic model variables that are relevant for the heat budget? • What can we learn about the North Sea system from the observed responses to applied temperature perturbations?
The exchanges of mass and heat among the North Sea, its adjacent seas, and the atmosphere are complex and influenced by processes on various time scales. The water temperature in the upper ocean is an important component in this system and leads to the main objective of our study, that is, to investigate and understand the impact of SST assimilation on physical processes that drive volume and heat transport.
In the present study, the 3-D numerical circulation model Nucleus for European Modelling of the Ocean (NEMO) coupled with the wave model WAM (Staneva et al., 2018) is employed and used for assimilation of Ocean and Sea Ice Satellite Application Facility (OSI SAF) SST satellite measurements (a detailed description is provided in the next section). Wave coupling is important for the temperature evolution of the North Sea because the inclusion of ocean-wave feedback modifies the stresses at the water surface and improves the simulation of the turbulent mixing layer thickness (Lewis et al., 2019). The DA scheme applied herein is based on the 3-D variational DA (3DVAR) analysis technique, which has frequently been employed in previous ocean studies (Dobricic et al., 2005;Liu et al., 2009). The complete 3-D temperature field is analyzed based on SST observations and assumptions about the vertical and horizontal structures of model errors.
The corresponding model error covariance matrix is not constant, but varies over time and space as a function of the mixed layer thickness. Regarding the heat budget of the North Sea, the analysis is performed as follows: • The total vertical heat fluxes, including latent and sensible heat fluxes, as well as short and long wave radiation components between the atmosphere and ocean are considered • The lateral advection of heat through open boundaries towards the Atlantic and Baltic Sea is analyzed along transects (see Figure 1) that were examined in previous studies (e.g., Hjøllo et al., 2009) • The tidal and wind-driven components of heat advection as well as those resulting from their interactions are considered separately In the assimilation of satellite SST observations, the response of the system to SST changes is analyzed such that the consequences of variations in the temperature, surface elevation, and current on the heat fluxes are decomposed into separate terms, accounting for different coupling mechanisms. Hence, the lateral advection of heat across the open boundaries of the North Sea is attributed to components related to Eulerian transport, tidal and wind-driven currents, and higher-order nonlinear interactions.
The remainder of this paper is organized as follows: Section 2 includes a description of the coupled NEMO-WAM model, the observational data used and the 3DVAR scheme; furthermore, this section also explains the experimental setup and methods used for analyzing the volume and heat budgets of the North Sea. The modeling results with and without DA are presented in Section 3, and the model results are further analyzed and compared with that of existing studies in Section 4. The physical processes that induce the advective transport of volume and heat through the open boundaries of the North Sea and the impacts of SST assimilation are discussed in Section 4 as well. The main conclusions are given in Section 5.

Model System
The Geesthacht Coupled cOAstal model SysTem (GCOAST) is built upon a flexible and comprehensive coupled model system that integrates the most important key components of regional and coastal systems, enabling the inclusion of information from observations (Ho-Hagemann et al., 2020;Lewis et al., 2019;Staneva et al., 2019). GCOAST encompasses (i) atmosphere-ocean-wave interactions, (ii) dynamics in and fluxes across the land-sea transition zone, and (iii) coupling of the marine hydrosphere and biosphere. The GCOAST integrates a circulation model (NEMO), wave model ( (Figure 1). The entire model domain is divided into 900 subdomains (approximately 100 × 100 km in size), with each run in parallel. The vertical grid is the NEMO standard σ-z* hybrid grid with 50 levels. The model time step adopts the time splitting methodology with a baroclinic time step of 100 s. Vertical turbulent viscosities/diffusivities are calculated using the Generic Length Scale (GLS) turbulence model (Umlauf & Burchard, 2003) with the "k-ϵ" closure scheme and the second-order algebraic model of Canuto et al. (2001).
The WAM model used here is based on the statistical description of wave conditions in the directional frequency spectrum at each active model grid point within a certain model area. The wave action conservation equation, which is complemented with a suitable description of the relevant physical processes, is used to follow the evolution of each wave spectral component. A detailed description is given by the WAMDI group (Günther et al., 1992;Komen et al., 1994;The Wamdi Group, 1988) and Janssen (2008). WAM Cycle 4.7, which is used for the GCOAST wave hindcasts, is an update of the former WAM Cycle 4. The basic physics and numerical code are kept the same in the new release. The source function integration scheme of Hersbach and Janssen (1999) and the model updates by Bidlot et al. (2007) are incorporated. Similar to the circulation model, the entire model domain is divided into 270 subdomains and run in parallel. The spatial resolution, regional coverage and meteorological forcing are the same as those in NEMO. The model and its performance for the study area are described by Staneva et al. (2016) in more details. Ocean waves influence circulation through a number of processes, including turbulence due to breaking and non-breaking waves, the transfer of momentum from breaking waves to currents in deep and shallow water, wave interactions with planetary and local vorticity, and Langmuir turbulence. The NEMO ocean model has been modified to take into account the following wave effects as described by Staneva et al. (2017); Staneva et al. (2019), Alari et al. (2016), and Wu et al. (2019): (1) Stokes-Coriolis forcing; (2) sea state-dependent momentum flux; and (3) sea state-dependent energy flux. Details of the NEMO-WAM model, boundary forcing and parameter settings are further explained in Appendix A.

Satellite SST Data
In this study, the modeled seawater temperature is analyzed with OSI SAF SST data, which are produced by the European Organization for the Exploitation of Meteorological Satellites (EUMETSAT). In this study, the North Atlantic Regional (NAR) SST is used, which consists of Metop/Advanced Very High Resolution Radiometer (AVHRR) and Suomi National Polar-orbiting Partnership (SNPP)/Visible Infrared Imaging Radiometer Suite (VIIRS) derived subskin SSTs over the North Atlantic and European seas at a 2 km resolution. In the study period (2017), these data are recorded twice a day, that is, at 10:00 and 20:00 UTC, with a different number of available points and at different locations due to varying cloud conditions. The SST products are compared with Match up Data Bases (MDB) of in situ (buoys) measurements, and they are classified into six levels at the pixel level: 0: unprocessed, 1: bad/cloudy, 2: worst, 3: low, 4: acceptable, and 5: excellent. For qualitative use, only values of levels 4 and 5 are employed in this study. More information on OSI SAF products can be found at http://osi-saf.org.

In Situ Observations
Another important source of information for the present study is fixed measurements from stations in the German Bight operated by the German Federal Maritime and Hydrographic Agency (Bundesamt für Seeschifffahrt und Hydrographie, BSH). Their Marine Environmental Monitoring Network in the North Sea and Baltic Sea (MARNET) consists of six automatic oceanographic stations in the North Sea, five of which are currently operating. Most stations measure temperature, salinity, oxygen, sea level, air temperature, wind direction, wind speed, and air pressure. In the present study, the two MARNET stations Nordseeboje III (NSB III, 54°41′ N and 6°47′ E) and FINO-1 (54°00.892′ N and 6°35.258′ E) (see Figure 1 for the locations) are used for validation purposes. All data are collected continuously on the platform and transferred hourly to the shore for further processing.

Assimilation Scheme
The assimilation technique used in this study follows the 3DVAR approach described in Lorenc (1997). The 3DVAR analysis scheme is based on the minimization of the cost function J defined as follows: where x f is the prior state, y is the observation vector, B is the model error covariance matrix, G is the observation error covariance matrix, and H is the observation operator. The most critical component of this approach is the definition of matrix B, which determines how the observation information is transferred to model regions where observations are not available. The most common approach is to prescribe correlation lengths for different dimensions in the model (e.g., horizontal and vertical). This matrix can also be used to consider physical relationships between different model variables to make the analysis dynamically consistent. As is usually the case in the existing literature, matrix G is assumed to be diagonal; that is, the observation errors are assumed to be spatially uncorrelated. An observation error standard deviation of 0.6 K was used, which is consistent with previous studies (e.g., Grayek et al., 2015). Because of the large state vector dimension, the explicit storage or the inversion of the a priori error covariance matrix B is prohibitively expensive. Therefore, one common approach (Lorenc, 1997) is to define a transformed state x as with the matrix C given as follows: The matrix C −1 can be thought of as an operator that removes correlations, because x has a diagonal covariance matrix. With these definitions, the cost function J becomes Here, the role of operator C is to add correlations. The gradient of J is given by or written as a column vector which can be re-formulated as The minimum of the cost function can be found using the condition ∇J(x) = 0, which is solved using a conjugate gradient (CG) method.
The CG method is an iterative technique that requires multiple applications of the matrix A (Equation 8) to given state vectors. Because of the high state dimension, these operations are computationally demanding. The 3DVAR implementation of the Geesthacht Assimilation System (GALATON3DVAR) applied in this study therefore uses the parallelized implementation of the circulation model. The filter operations are distributed over different processes assigned to the various model subdomains. The required information exchange between subdomains is implemented using the Message Passing Interface (MPI). GALATON3D-VAR uses pre-computed index tables and is applicable for both regular and irregular model grids.
In the standard 3DVAR scheme, the error covariance matrix B is constant over time. In this study we make the vertical error correlation length, which is included in the formulation of B, dependent on the mixing layer thickness. This approach is based on the assumption that the mixed layer thickness of the free model is within a reasonable error margin. In that case, by definition, the complete mixed layer temperature will be affected by the same error as the SST. Thus, the vertical model error correlation length varies both in time and space with the evolution of turbulent mixing. In technical terms, we have implemented the vertical error correlations using a filter C of the form where σ is the mixed layer thickness, x i represents the temperature at the ith layer and z i is the respective distance from the surface. If we start with a white noise process x with variance   2 0 x , the filtered signal has a covariance that is, the error correlations with the surface temperatures are decaying exponentially with increasing depth. For the present study a standard deviation of SST model errors of   2 0 x = 0.5 K was assumed. The assumptions about vertical error correlations represented by Equation 11 has to be seen as a compromise to follow two objectives, which are contradicting to a certain extend.
(1) The analysis should avoid unrealistic corrections of the stratification; (2) The model error correlations with the surface errors should decay with depth. The functional form used to describe error decorrelation allows to define a vertical scale of information propagation from the surface (e-folding distance), but one has to accept that the stratification of the analysis is not in all cases fully consistent with the Free Run (e.g., by introducing small density gradients in the mixed layer). To further optimize the correlation model Equation 11 would require a much deeper analysis of model error characteristics, which is beyond the scope of the present analysis.
Furthermore, the horizontal correlation length of model errors required for the formulation of matrix B is assumed to be constant over both time and space with a value of 5 km. This value represents a compromise between retaining small-scale structures introduced by the satellite measurements into the system and simultaneously extrapolating observations to fill data gaps. Of particular concern in this context were eddy structures (e.g., along the shelf edge), which are known to play a role in horizontal heat transports.
It is important to emphasize that the mixed layer thickness is significantly dependent on ocean waves, which impact both the air/sea momentum fluxes and the turbulent kinetic energy in the upper ocean layers. For this reason, the coupling of the circulation model and the ocean wave model is an important factor in the applied DA method.
For a multivariate 3D-Var assimilation, sea surface height (SSH)-seawater temperature balance is also considered, following the approach of NEMOVAR system developed by ECMWF (Mogensen et al., 2012). The baroclinic SSH component δη is linearly balanced with changed seawater temperature δT, and estimated by computing the dynamic height at the water surface relative to the mixed water depth: where ρ r is a reference density, and Here, α is the thermal expansion coefficient and is estimated through a direct tangent-linearization of the equation of state (EOS) (Jackett & McDougall, 1995). We use the EOS formulation as given in and the density of sea water (at 1 atm) is given by Millero and Poisso (1981) In this equation, S is the salinity of sea water in ppt (parts per thousand by volume). The reference density ρ r , the coefficients A, B, and C are also functions of temperature. For the equations and further details regarding balance operations, the readers are referred to Appendix C.

Experimental Setup
To initialize the DA run with a relatively balanced state, a spin-up run of more than 3 years (from October 15, 2013 to January 1, 2017) is performed with the coupled GCOAST system. Subsequently, the coupled model continuously runs for 1 year without assimilating OSI SAF SST into the system. This case is denoted as the "Free Run" case. Then, the model is restarted from January 1, 2017 with the same settings, but the OSI SAF SST data are assimilated; this run is denoted as the "DA Run." In the DA Run, only the satellite data at 10:00 each day are used for the analysis according to the 24-h assimilation interval. The SST assimilation is performed at noon. The analysis time has about 2-h time separation with respect to the observations and potential temperatures changes within that time window are considered so small that the results are very close to a First Guess at Appropriate Time (FGAT) implementation. The primary reason for the analysis at 12:00 is that this is the main synoptic hour closest to the observations. The main synoptic hours (0, 6, 12, 18 UTC) are the standard times used in operational forecast systems to carry out an analysis and to launch a new forecast. A typical amplitude scale for the diurnal SST cycle is 1°C with a maximum occurring in the afternoon around 15:00. This means that the expected temperature variation between 10:00 and 12:00 is below 0.2°C, which is significantly smaller than the scale assumed for the observation errors (compare Section 2.3). As temperature is usually increasing at this time of the day, this deviation is systematic, however one also has to put this into perspective with the total temperature variations covered in this study. The analyzed annual cycle shows typical temperature variations of the order of 10°C in the considered region, and the error introduced by the time mismatch is almost two orders of magnitude smaller than that.
Furthermore, "the DA Run" experiment is repeated but with NEMO uncoupled with WAM, and the same OSI SAF SST data set is applied for the assimilation. This configuration is used only for a comparison with the coupled model run to inspect the impact of coupling on the SST simulation. Hence, in the following sections, the DA Run refers to the coupled assimilation run unless otherwise specified as being "uncoupled." The set-up specified for each experiment is summarized in Table 1.
To determine the water and heat exchanges between the North Sea and the adjacent seas, five transects are selected through the western (Dover Strait), northern (Fair Isle, Shetland Shelf, and Norwegian Channel) and eastern (Danish Strait) open boundaries of the North Sea ( Figure 1); these transects are equivalent to sections N13, N3, N1 N2, and N22 of the North West European Shelf Operational Oceanographic System (NOOS), respectively (NOOS, 2010). The modeled sea surface elevation, velocities and seawater temperature are output at 1, 3, and 6 h intervals, respectively. To compute the transport of volume and heat across these transects, the 3-hourly velocities and 6-hourly seawater temperature are both interpolated to hourly resolution with a "spline interpolation" method.

Volume and Heat Transports
To analyze the volume budget for the North Sea, the water transport through the five transects illustrated in Figure 1 is considered. The water transport is computed as where A is the area of the 2D transect plane and u is the current component perpendicular to the plane. Regarding the heat budget, both lateral heat transport and air-sea heat exchange should be considered, with the former due to the transport of water with different temperatures and the latter caused by four main processes: heat fluxes due to shortwave and longwave radiation and latent and sensible heat fluxes. Note that the present study focuses on the total heat budget over the North Sea and thus does not distinguish individual processes in the atmosphere-ocean heat flux exchange. The advective heat transport is given by and the incremental area is expressed as dA = dh dz, where dh and dz are the grid sizes in the horizontal and vertical directions, respectively. Note that due to the σ − z* hybrid grid used in the model, dz depends on the water elevation and local depth. The constant c p = 4.19 × 10 3 J kg −1 K −1 is the heat capacity, ρ = 1,026 kg m −3 CHEN ET AL. is the reference seawater density, and T K is the temperature in K. Note that positive flux values refer to mass and heat flow into the North Sea; that is, positive is northeastward along the transect Dover Strait transect, southeastward along the Fair Isle transect, southward along the Shetland Shelf, and Norwegian Channel transects, and westward at the Danish Strait transect. Hence, summing all the positive and negative transports of grid cells separately over the study period yields the inflow and outflow through each transect. Then, the net transport is obtained by adding up the total inflow and outflow. This approach is adopted to avoid large numerical rounding errors that may occur in a straightforward sequential summation of position and negative numbers.
Heat transport is affected by changes in the water level, current speed and temperature, and all of these effects and their coupling will be analyzed in the following. First, the velocity normal to the transect, the water layer thickness and the temperature in the DA Run are written as u a = u f + u d , dA a = dA f + dA d , and T a = T f + T d , respectively, with the subscript d denoting the differences in variables between the Free Run (f) and DA Run (a). Then, following previous studies (Dieterich et al., 2019;Graham et al., 2016;Hjøllo et al., 2009), the heat transport is computed with respect to a reference temperature T rK . In our analysis, we adopted the approach taken in Dieterich et al. (2019) and applied a reference temperature of 6°C (i.e., T rK = 279.15 K). Thus, Equation 16 is decomposed for the DA Run as follows: where V 1 is the Free Run volume transport, V 2 is the volume transport due to the modification of the velocity, V 3 is the volume transport due to the modification of the sea surface height, and V 4 is associated with the nonlinear interaction between the changed velocity and sea surface height. Likewise, the heat transport (Equation 17) for the DA Run reads As shown in Equation 19, the DA Run heat transport consists of Part I: The heat transport relative to TrK for the free run Part II: The change in heat transport due to the DA-induced temperature change; Part III: The change in heat transport due to the DA-induced velocity or area change; Part IV: The change in heat transport due to nonlinear interaction between the DA-induced temperature changes, velocity changes and area changes; Part V: The change of reference heat transport.
The decomposition of H a q distinguishes the direct changes (Part II) and the secondary changes (Part III and IV) in advective heat transport due to assimilation of temperature. Moreover, it is worth to note the difference between the heat transport relative to the reference temperature of the Free Run (P 2 ) and the change of the reference heat transport (because of changed velocity and area) in the DA Run (Part V). Thus, for a consistent quantification of the heat transport relative to T rK for the DA Run, P 4 + P 10 , P 5 + P 11 and P 6 + P 12 obtains the influence of DA due to velocity changes, area changes and their nonlinear interactions, respectively.

Improved Consistency of the Modeled Seawater Temperature
Before evaluating the transport and budget in the North Sea, the general DA performance is first assessed in the study period for the GCOAST domain. Figures 2a and 2e show the coverage of the OSI SAF SST in the warming season (April to September) and the cooling season (January to March and October to December). Although variations are observed in both space and time (due to cloud cover), the OSI SAF data for the North Sea are available for 30% ∼ 40% of the time during the warming season and for 10% ∼ 20% of the time during the cooling season. Figures 2b, 2c, 2f and 2g are the Free Run and DA Run (analysis) on 1 July and 31 December 2017, which are in the middle and at the end of the assimilation period, respectively. The model updates introduced by DA are kept and accumulated over time. Temperature differences between the model and the observations are small in July when compared to that in December. In the North Sea, the SST is increased by 0.5°C-1°C due to DA (Figure 2h). In the Baltic Sea, the Free Run SST is lower than the DA Run SST in July (Figure 2d), while the former is higher in December (Figure 2h). Furthermore, it is worth noting that the large spatial scale SST features in the Free Run (Figures 2b and 2f) are kept in the DA Run (Figures 2c and 2g); moreover, smaller-scale features are introduced by the analysis (e.g., comparing the SST structures in Figures 2f and 2g at 60 to 62°N, −5 to 5°E on 31 December). Note that the OSI SAF SST is shown at 10:00 while the model data are shown at 12:00.
The numbers of available OSI SAF data points for the model grid used for DA over the full model domain are shown in Figure 3a. Due to cloud cover, the number of valid observation data points roughly varies from 1 × 10 4 to 7 × 10 4 . In general, the number of valid observations is smaller in the cooling season than in the warming season. Figure 3 also compares the temporal evolution of RMSEs for the analyzed states of the model with and without DA, and the results show that the Free Run has an RMSE of 0.7°C ∼ 1.8°C, whereas after DA, the RMSE is reduced to 0.3°C ∼ 1.2°C. Furthermore, the DA scheme keeps the analyzed SST close to the observations even when the Free Run shows strong departures. For example, in the early stage of the study period (January to June), the RMSE of the DA Runs is close to that of the Free Run (difference is approximately 0.1°C ∼ 0.3°C), while at the end of 2017 (October to December), the difference between CHEN ET AL. the Free Run and the DA Runs RMSE becomes 0.5°C ∼ 0.6°C. Compared with the model runs with and without DA, the RMSE difference of SST between the coupled and uncoupled DA Runs is small. Figure 3 shows that the analysis is able to bring the model fields closer to the observed data. However, the objective is not to achieve perfect agreement because the observations are known to contain errors. The remaining RMS difference between the analysis and observations is consistent with the order of the assumed observation error, thus indicating that the assumed measurement errors are reasonable. Figure 4 compares the modeled (coupled DA) temporal variation of temperature with that of the independent MARNET data acquired at NSB III and FINO-1 (see Figure 1 for their locations) at different water depths. The OSI SAF SST data over the study period are plotted at these two locations together with the near-surface temperature data at a depth of 3 m. Note that the OSI SAF SST is unavailable when the station is covered by clouds or the data quality is low. The comparison at this water layer verifies the consistency between the OSI SAF data and in situ data. At the NSB III station (Figures 4a-4c), the Free Run has already shown a good capability to simulate the temperature evolution in 2017 despite an underestimation of approximately 1°C-2°C in the middle water layers between June and July (Figures 4d-4f). It is worth to note that the impact of DA is not always positive. For instance, at NSB III, a deterioration is observed between July and October near the bottom (Figures 4c and 4f). At FINO-1, the DA leads to improvements in the performance of seawater temperature simulation (Figures 4g-4i). During the first half of the study year (January to July), the underestimated temperature (approximately 1°C) between the Free Run and the MARNET data is corrected in the DA Run. However, when the water column starts to stratify (in the beginning of summer) the temperature at 15 m of the DA Run shows similar behavior as one of the Free Run (Figures 4b and 4e). Due to the relatively shallow water depth and mixing of the water column, the DA Run also shows representative temperature evolution near the sea bed (Figures 4j-4l). The temperature difference between the CTD measurements at FINO-1 at the last quarter (October-December) and both runs shows similar behavior. CHEN ET AL.

Volume Budget
The volume transport through each of the five selected transects in the North Sea (see Figure 1) is computed over the study period (see Equation 16) and illustrated in Figure 5a. The transects are selected following previous studies (NOOS, 2010). Net transport into the North Sea is observed through the Dover Strait, Fair Isle, Shetland Shelf, and Danish Strait, while net outward transport occurs through the Norwegian Channel. In the Free Run, the net transport through the Danish Strait (approximately 0.04 Sv, 1 Sv ≡ 10 6 m 3 s −1 ) CHEN ET AL.   is the smallest among the inward transects, followed by the Dover Strait (0.088 Sv). In contrast, at Fair Isle and Shetland Shelf, the net transport is around 0.48 Sv, accounting for 44% of the net water transport into the North Sea, respectively. Regarding the outward volume transport, 1.27 Sv exits through the Norwegian Channel. In the DA Run, the inward volume transports are 0.076 Sv (7.4%) through the Dover Strait, 0.45 Sv (43.5%) at Fair Isle and 0.47 Sv (45.2%) along the Shetland Shelf and 0.04 Sv (4%) the Danish Strait, while the outward transport is 1.19 Sv through the Norwegian Channel. The total volume transport into the North Sea is approximately 15% (Free Run) and 13% (DA Run) lower than the total outward transport, respectively. Compared with the volume exchange between the North Sea and the Atlantic, river runoff is negligible. The total river discharge into the North Sea is approximately 300 cubic km/year (Quante & Colijn, 2016), which is only approximately 0.01 Sv. The modeled volume budget (as well as the heat budget presented in the next section) in the present work are further compared with those of existing studies in Section 4. Figure 5b compares the volume transport composition of the DA Run along each transect, demonstrating that the volume transport differences between the Free Run and DA Run are associate with changes in the current field (V 2 in Equation 18), whereas the other contributions are negligible. Interestingly, the sign of the volume transport due to modified velocities is, except for the Danish Strait, opposite to that of Free Run (V 1 in Equation 18), which implies that DA is mostly decreasing the volume transport. The potential physical processes responsible for this are further explored and discussed in Section 4.

Heat Budget
Neglecting the small heat flux associated with river runoff, the heat content in the North Sea depends mainly on two processes: lateral heat transport through the open boundaries (the five selected transects illustrated in Figure 1) and heat exchange between the air and sea at the water surface. Figure 5c shows the heat transport through each of the selected transects. Since the advective heat transport is largely determined by the volume transport, the heat transport distribution through the open boundaries of the North Sea is similar to the volume transport distribution; that is, the North Sea gains heat through the Dover Strait, Fair Isle, Shetland Shelf, and the Danish Strait, loses heat through the Norwegian Channel. However, due to differences in the local temperature features of seawater, the heat transport ratios through these transects are different from the volume transport ratios. Despite the low volume transport through the Dover Strait, the high seawater temperature results in a large amount of heat into the North Sea (2.67 TW in the Free Run, which accounts for 14.55% of the total transport, and 2.59 TW in the DA Run, i.e., 14.48%). Similarly, through the Danish Strait, more cold water transported inward the North Sea compensated warmer sea water outflow to the Baltic, thus yielding a negligible heat transport through the Danish Strait. In total, the North Sea gains 18.32 TW (17.78 TW) of heat and loses 16.32 TW (15.04 TW) in the Free Run (DA Run), with a net of 2.02 TW (2.74 TW), which yields net heat gain/heat loss ratios of 1.12 and 1.18 for the Free Run and DA Run, respectively. In other words, with DA, more heat enters the North Sea via lateral advection, during the considered period. Further comparison of the current work and existing studies are discussed in Section 4.1 Decomposition of the heat transport for the DA Run (Figure 5d) reveals that the heat transport relative to T rK for the Free Run is the main component (Equation 19 Part I, P 1 + P 2 ). Note that the choice of a reference temperature will change P 2 and thus the relative composition, which is further discussed in the next section. The direct change in heat transport between the DA Run and the Free Run due to temperature difference (Equation 19 Part II, P 3 ) is only significant in the Danish Strait. Other terms that contribute to the DA-induced heat transport are those result from the changes in volume transport due to DA (Equation 19), in which P 4 + P 10 , the heat transports due to the DA-induced velocity change, is the major contribution. As discussed in the previous section, this term is correlated with the secondary effect of temperature assimilation on hydrodynamics. The nonlinear interactions between the DA-induced temperature changes, velocity changes, and area changes (Equation 19 Part IV) are negligible in contributing to the changes in heat transport.
With regard to air-sea heat exchanges, Figure 6 shows the annual mean heat fluxes over the North Sea. Evidently, the assimilation of SST data does not change the spatial pattern of air-sea heat fluxes. In both the Free Run (Figure 6a) and the DA Run (Figure 6b), the North Sea gains heat from the atmosphere in the middle of the domain close to the British coast and loses heat along the northern and southern boundaries.
The Free Run yields a net heat gain for the central North Sea with a maximum value of 30 W m −2 and a net heat loss for the area close to the northern boundary, with a minimum value of −50 W m −2 . In the DA Run, the maximum net heat gain is approximately 10 W m −2 , while the maximum net heat loss is approximately −60 W m −2 . As shown in Figure 6c, the assimilation of SST data results mainly in a reduction in the heat gain from the atmosphere in the middle of the North Sea and the Norwegian Channel. For the entire North Sea domain, the total annual mean heat flux from the sea to the air has increased, with the value change from 7.97 to 10.50 TW (i.e., a constant air-sea heat flux change of −4.8 W m −2 ). Figure 7 shows the annual accumulation of net heat transport over the North Sea. The heat transport between the air and sea exhibits a strong seasonal cycle (approximately 3.5 × 10 20 J), with warming of the North Sea from middle March to September and cooling during the remainder of the period. The heat energy accumulated through ocean-atmosphere fluxes is negative at the end of the year, which means that there is a net loss of heat from the North Sea to the atmosphere over the year 2017. Advective heat transport shows a continuous heat gain trend from the North Atlantic Ocean to the North Sea. In the Free Run, the North Sea loses approximately 2.5 × 10 20 J via the air-sea interface, while it gains 0.62 × 10 20 J from lateral water transport. In the coupled DA Run, the North Sea loses approximately 3.3 × 10 20 J via the air-sea interface and gains 0.89 × 10 20 J from lateral water transport. In other words, both heat transport processes are enhanced by DA, albeit with a different sign. This finding indicates that, with DA the North Sea gains more heat through lateral advection (0.27 × 10 20 J) and loses extra heat due to air-sea heat exchange (0.8 × 10 20 J). CHEN ET AL. Thus, the net heat gain by lateral advection warms the North Sea, while the increased SST by assimilation increases the heat gradient between the sea and the atmosphere, which further enhances the heat flux between them. A comparison of the two experiments shows that the amount of net heat loss due to air-sea heat exchange compensated by advective transport remains more or less the same (with the ratio of 25%-27%). Wave coupling also contributes to enhancing the net heat gain in the North Sea via advective transport (by modifying the momentum budget at the ocean surface [Lewis et al., 2019]). It is worth to stress that despite the impacts of wave-circulation coupling are minor for the SST (see Figure 3), it is important when considering the water temperature distribution in the vertical direction. As it is shown in Figure 7, during the warming period, the advective heat transport in uncoupled model run (solid black curve) is lower when comparing with the coupled run (solid red curve). Moreover, the uncoupled model run also lowers the air-sea exchange, especially between August and November.

Transport Estimates
To assess the transports computed in the current study, earlier investigations are reviewed regarding the average rates of volume and heat exchanges between the North Sea and its adjacent seas. The net inflow through the Dover Strait to the North Sea obtained in the present study was approximately 0.076 ∼ 0.088 Sv, which is close to the field measurements (0.09 Sv) reported by Prandle et al. (1996), who estimated these values based on high-frequency (HF) radar and bottom-mounted acoustic Doppler current profiler (ADCP) profiles spanning 1 year. At Fair Isle, the net inflow of 0.48 Sv in the Free Run (0.45 Sv in the DA Run) is higher than the observations (0.36 Sv) reported by Otto et al. (1990) but close to the value of 0.49 Sv obtained by taking the annual mean of a model study (Winther & Johannessen, 2006). A consistent value is also found for the Shetland Shelf, where the observed net inflow is approximately 0.6 Sv (Otto et al., 1990). At the Norwegian Channel, Otto et al. (1990) observed an inflow flux of approximately 0.7-1.1 Sv and an outflow flux of 1.8 Sv, thus yielding a net outward flux of 0.7 ∼ 1.1 Sv from the North Sea. Based on the model study of Winther and Johannessen (2006), Hjøllo et al. (2009) computed the annual mean volume transport through the Norwegian Channel and found a value of 1.23 Sv (2.33 Sv) into (out of) the North Sea. Hjøllo et al. (2009) Omstedt et al. (2004). There are two main water masses: the inflowing Norwegian Coast Current (NCC) on the northern side of the transect and the outflowing North Sea Jutland Current (JC) on the southern side, and are roughly balanced in exchange (Danielssen et al., 1997;Hjøllo et al., 2009;Winther & Johannessen, 2006). Although the modeled water transport in time series differs in the details, the findings illustrate similar patterns when further compared with the BSH model, which was applied to the same transects (Dick et al., 2001).
Although the total heat content in the North Sea is independent from the reference temperature, one should note that the choice of the reference temperature (T r ) is rather arbitrary. The purpose of using T r , related to data assimilation, is to allow an easy comparison of the relative changes. Thus, the amount of heat transported by lateral advection of the water mass for the temperature higher than T r is considered "active" heat content. Another concern in the choice of T r is the consistency with previous studies (e.g., T r = 0°C in Hjøllo et al., 2009, T r = −1.95°C in Graham et al., 2016, andT r = 0°C in Dieterich et al., 2019). In the following, the sensitivity of the total heat budget to the reference temperature (T r ) is investigated by computing heat transport compartments (Equation 19) with different given T r . Figure 8a shows the net heat transport through the transects for T r = 0°C. As expected, the net heat transport through each transect is much larger than that in the T r = 6°C case (see Figure 5a), since P 2 in which T r greater than 0°C reduces the transport term P1. In this case, the total advective heat transport is −2.9 TW in the Free Run (−1.7 TW in the DA Run).
Since the advective heat transport into the North Sea exhibits a strong annual variation (Hjøllo et al., 2009), it is difficult to compare the exact values of the present work with those of early studies. However, the value range is comparable in magnitude to the results of our study. For example, Hjøllo et al. (2009), who also used T r = 0°C, computed the heat transport through the western boundary (0°E, 49 ∼ 50.5°N), the northern boundary (2°W ∼ 5°E, 59.2°N), and the eastern boundary (8.1 ∼ 8.6°E, 57.1 ∼ 58.1°N), and the results showed that the monthly mean net transport in 1985 ∼ 2007 through these boundaries varied between 0 and 10 TW, −15 and 8 TW, and −5 and 8 TW, respectively. The total North Sea advective inflow heat reached a mean of 2.6 TW during this period. Using the monthly means of net advective heat flux, we estimated the annual mean advective heat flux of about 0.5-7 TW during the study period of Hjøllo et al. (2009).
The choice of T r affects the value of P 2 in Equation 19, resulting different "active" advection of heat transport for the Free Run. Moreover, changing T r further changes the reference heat transport Part V of Equation 19). Figure 8b illustrates the sensitivity of DA transport changes to the choice of T r . It is shown that the amplitude of reference heat transport (green dashed line) increases when larger T r is applied. With the fact that the reaction in heat transport due to the DA-caused temperature change (cyan dashed line) is independent of T r , the DA leads to an increase of the "active" heat transport for small reference temperatures but a decrease for high reference temperatures. Figure 7a shows that the lowest heat content due to air-sea heat exchange occurs in late March, while the highest heat content exists in early September for both the Free Run and the DA Run. This finding implies that the DA mainly changes the amplitude of air-sea heat exchange but hardly has an effect on the temporal evolution pattern itself. For both the Free Run and the DA Run, the seasonal variation of the air-sea heat exchange is approximately 6.5 × 10 20 J (the difference between the net heat transport maximum in September and the minimum in late March). This finding is consistent with the 20-year model study result of Hjøllo et al. (2009), who defined the seasonal variation of the air-to-sea heat transport as a heat gain of 5 ∼ 6 × 10 20 J. Meanwhile, the net heat gain and heat loss ratio in our 2017 case (1.12 for the Free Run and 1.18 for the DA Run) is close to the calculation results from their study, which varied between 0.92 and 1.23.

Mixing Layer Thickness
As explained in Section 2.3, the thickness of the mixing layer is used to define the correlation length for the temperature assimilation. In this sense, the mixing layer thickness provides a depth scale over which the observation information is distributed. A well-resolved mixing layer depth is necessary to eliminate systematic errors in ocean temperature simulations. Earlier studies have identified the need for wave coupling to improve the performance for the annual variation of the mixing layer depth, in which the impact of coupling on temperatures could be diminished by DA (Lewis et al., 2019). This finding also motivates the implementation of the coupled NEMO-WAM model in the present study. The results show that wave CHEN ET AL.  coupling has an impact (e.g., Figures 3 and 7) for water temperature profiles in depth, particularly when the North Sea is rapidly warming up and cooling down. Nevertheless, the circulation-wave interaction is rather more complex than simply a weakening or enhancing of the mixing layer depth at a certain time or location. Details regarding the physical processes included in the coupling and their contribution to the ocean circulation are further described in Appendix B. However, the focus of this study is on the DA related aspects of estimating the heat budget; thus, we will perform a more detailed analysis of the coupling mechanisms in a separate follow-up investigation. Figure 4 shows that at NSB-III, the SST DA can correct the water temperatures in the middle water layers (15 m). However, DA hardly improves the temperatures in the deep water from May to July. This finding is consistent with the annual variation of the mixing layer thickness, which is only 10 m during this period and reaches the bottom during the remainder of the year (not shown). Under well-mixed conditions, such as at FINO-1 and along the Dover Strait and Fair Isle transects, the assimilation of SST data can improve the temperature over the entire water column.
To further discuss the behavior of the mixing layer depth and the associated temperature profile corrections in partially stratified water columns, Figure 9 illustrates the annual variation of the seawater temperature, including Free/DA Run difference, and mixing layer thickness at the Shetland Shelf and Norwegian Channel. Because the difference in layer thickness between the Free Run and the DA run is small, only one experiment (Free Run) mixing layer depth is plotted. At these two locations, the water column is stratified in the summertime and becomes increasingly mixed as the water temperature in the upper layer cools. The thickness of the temperature correction agrees well with the mixing layer depth distribution. In the cooling period (September to November), the largest changes are visible below the mixed layer depth. Moreover, the water layer "memorizes" the corrections of the prior state. Such a memory could last for months if the hydrodynamic conditions are sufficiently stable. For example, at the Shetland Shelf (Figure 9b) in early May, the temperature correction reaches 100 m above the bottom, which is consistent with the mixing layer thickness. From late May until September, the stratification is rather stable and the mixing layer thickness remains only 25 m. Within this layer, the water temperature is reduced (−0.5°C) due to the DA, whereas in the layers below down to 100 m above the bottom, the temperature corrections (approximately + 0.5°C) performed in early May are still visible.
DA also cools the surface temperature at the Norwegian Channel (within the top 25 m between June and late August). However, at water depths below the mixing layer depth, the temperature is corrected in the opposite direction. From late July to early August, the temperature of the water column between 40 and CHEN ET AL. In Equation 23, two additional terms occur, which are related to transports caused by tidal pumping and wind pumping. The former is due to the correlation between seawater temperature fluctuations and tides, while the latter is due to the correlation between seawater temperature fluctuations and wind.
As shown in Figure 10, along all the transects throughout the North Sea, the main contribution to volume exchange is Eulerian transport. Stokes transport, which results from the correlation between strong tidal currents and tidal waves distorted by bottom friction, also plays an important role at the Dover Strait (contributing approximately 0.04 Sv). Along all the remaining transects, Stokes transport, wind and wind-tide interactions are weak. The results reveal a considerable difference between the Free Run and DA Run at the Norwegian Channel due to the strong reduction in Eulerian transport. Moreover, with regard to heat advection ( Figure 10, lower panel), Stokes transport is as significant as Eulerian transport through the Dover Strait. At the northern boundaries, (Fair Isle, Shetland Shelf, and Norwegian channel), the net transport into the North Sea due to the annual mean current is an order of magnitude larger than that due to the other processes. At the Norwegian Channel, tidal pumping and wind pumping of heat are also important. In the Danish Strait, tidal pumping and wind pumping are equally important but play opposite roles, while the Eulerian transport is negligible. Similar to the results of the volume decomposition, at the Norwegian Channel, DA reduces heat transport by decreasing the annual mean current.
Because the atmospheric forcing is prescribed, changes in the annual mean wind stress are not considered. The main differences in Eulerian transport between the Free Run and DA Run are attributed to modifications, such as spatial variations in heating/cooling of the surface, which lead to changes in (a) pressure gradients and (b) stability of the water column and hence reduced/increased internal friction. Basically, an increased stability of the water column will reduce internal friction and vice versa. Nonetheless, the resulting current modifications are of course complicated, because they are also controlled by other factors like bathymetry.
CHEN ET AL.  Figure 11 shows the annual mean sea surface height (SSH) difference and density difference between the DA Run and the Free Run along the Norwegian Mouth transect (see Figure 1). A high correlation is observed between the two terms, with reduced water density inside the Norwegian Channel (at distances of less than 100 km, Figure 11b) and increased density near the shelf edge (120 ∼ 200 km, Figure 11b). Furthermore, as mentioned in the previous section, an inward current from the Atlantic to the North Sea exists along the western side of the Norwegian Channel transect; this current represents a branch of the eastward along-shelf current in the Atlantic. This finding implies that the observed modifications of the pressure gradient and thus the current field must be due to changes of internal friction caused by adjustments of the temperature stratification.
In Figure 12, the arrows with length and direction reflect the vertical mean current averaged over monthly periods (hence, the dominant tidal components are removed, e.g., M 2 ) for the different seasons of 2017. The contours indicate the relative difference in the dynamic energy between the DA Run and the Free Run per unit water mass that is,   Here, red (positive values) refer to the enhancement of the current due to the assimilation of SST data, while blue (negative values) indicates weakening of the current. When ɛ = 1, the water current is generated fully by the DA, whereas for ɛ = −1 DA removes water current locally. The case ɛ = 0 implies that DA does not change the local velocity field. Figure 12 clearly shows that the larger changes of the vertical mean current due to DA occur in area with deeper water depth. The largest dynamic changes occur inside the Norwegian Channel and outside of the shelf sea. Currents in the Atlantic are enhanced along the northern side of the North Sea shelf because of DA (denoted by the yellow arrow in Figure 12). As a result, the branch current along the western side of the Norwegian Channel increases persistently over the whole year. Such enhanced inflow compensates for the net outward flow of the North Sea.

Conclusions
The present study investigated the impact of the assimilation of satellite SST data on the simulation of the volume and heat budgets for the North Sea in a wave-circulation coupled system. This work follows that of Lewis et al. (2019), who showed that model coupling alone is not a sufficient strategy for improving all aspects of model performance. The 3DVAR scheme is implemented with the assumption that the model SST errors are strongly correlated with the temperature errors inside the mixing layer. This work demonstrates that the assimilation of OSI SAF SST data can improve the model analysis results as shown by the reduction CHEN ET AL.

10.1029/2020JC017059
19 of 25 in the RMSE from 0.7°C ∼ 1.8°C to 0.3°C ∼ 0.9°C. With DA, the overall modeled surface temperature in the North Sea increases by approximately 0.5°C ∼ 1.2°C. Further comparison of the analyzed model results with independent profile observations shows improvements in the seawater temperature at deeper water layers as well. At the two MARNET stations NSB III and FINO-1, the differences (1.0°C ∼ 1.5°C) between the observations and the Free Run at moderate water depths and near the bottom are mostly corrected by the assimilation, despite to a certain extend errors still remain.
The total lateral advective volume and heat transports of the DA Run are decomposed into terms that are induced by the individual responses of physical quantities and their interactions with the assimilation. These quantities include the current velocity, sea surface elevation and water temperature. With DA, both the inward and the outward volume transport and consequently the lateral heat transport of the North Sea are decreased. The main difference in the volume transport between the Free Run and the DA Run is associated with the current velocity changes induced by the assimilation of SST data. This term counteracts the Free Run volume transport and has the largest impact on the Norwegian Channel among all five open boundary transects crossing the North Sea.
The current field changes caused by SST assimilation further affect the lateral heat transport. The net heat gain and heat loss ratio of the North Sea due to advection has increased from 1.12 to 1.18 based on DA. As another main component of the heat budget of the North Sea, the air-sea heat exchange is enhanced due to SST assimilation. The results show that the heat flux from the sea to air increases by 31% (from 7.97 TW to CHEN ET AL. 10.1029/2020JC017059 20 of 25 10.5 TW) over the entire domain because of direct local temperature corrections by DA and indirect impacts from the non-local temperature correction and heat transport by lateral advection.
Further analysis of the advective heat exchange induced by individual physical processes reveals that Eulerian transport is the dominant mechanism along the northern (along the Fair Isle, Shetland Shelf, and Norwegian Channel transect) open boundary of the North Sea. In the west (the Dover Strait transect), Stokes transport is as important as Eulerian transport, while in the east (the Danish Strait transect), tidal pumping and wind pumping play an equal role. At the Norwegian Channel, heat transport due to tidal pumping and wind pumping have second-order contributions to the total heat exchange. For all the northern side open boundaries of the North Sea, heat transport induced by Stokes, the annual mean wind stress, and wind-tide interactions is rather small. SST assimilation decreases the horizontal density gradients from the North Sea to the Atlantic, causes changes in the velocity field within the Norwegian Channel and in turn reduces the Eulerian transport of water and heat from the North Sea to the Atlantic. As the analysis contains a balancing operator for pressure gradients introduced by the temperature modifications, the observed model response is dominated by modifications of the internal friction in the 3D circulation. The resulting changes in lateral heat fluxes are significant. Because of the observed coupling between lateral heat advection and DA corrections of the mixed layer. One limitation of the study is caused by missing observations of currents, in particular across the shelf edge to further constrain the heat budget. The analysis has also shown that the coupling between the circulation model and the wave model has a notable impact on the analyzed heat budget. This effect is expected more generally for SST assimilation systems, where the focus of the model update is on the mixed layer.
In summary, this study quantifies the changes induced by 3DVAR data assimilation in different components of the heat budget of the North Sea. With the growing debate regarding the impacts of climate change on the North Sea, accurate assessments of these factors will be of increasing importance in the future.

Appendix A: NEMO-WAM Setup
A side view of the σ-z* hybrid grid is sketched in Figure A1. The standard NEMO hybrid grid features tangential stretching below a depth of 200 m. The minimum water depth of the model is 8 m, and the maximum depth is 6,300 m. This grid configuration results in a minimum level thickness of 0.16 m at the surface and a maximum level thickness of 755 m at the bottom. The model uses a nonlinear free sea surface with a variable volume. The barotropic subcycle time step is adaptively determined during the run time using a maximum Courant number of 0.8. The solution of the barotropic subcycle is filtered using a box filter, and the hydrostatic pressure gradient is estimated using the σ-coordinate pressure Jacobian scheme.
The advection of momentum conserves both energy and enstrophy, and a bi-Laplacian horizontal diffusion operator with a coefficient of −2.8 × 10 −8 m 4 s −1 is applied. The free-slip lateral boundary condition is employed for momentum along the coastline. A spatially varying logarithmic layer is adopted for the bottom friction with a bottom roughness length of 1 × 10 −6 m for the Danish Straits and the Baltic Sea and 4 × 10 −3 m for the rest of the domain. Tracer advection uses a total variation diminishing (TVD) approach with the total variance decreasing scheme described in Zalesak (1979). Tracer advection is parameterized by the TVD scheme and a 2D varying Laplacian diffusion operator with a value of 0.25 m 2 s −1 for the region covering the Danish Straits and the Baltic Sea and a value of 50 m 2 s −1 for the rest of the domain. Lateral diffusion for the tracers is applied along with geopotential levels. In the vertical direction, the Craig and Banner surface wave mixing parameterization (Craig & Banner, 1994) is applied with a wave breaking TKE flux constant of 150, and the dissipation under stratification is limited using a Galperin limit of 0.07 (Galperin et al., 1988).
The lateral boundary forcing uses the NEMO "BDY" standard. The boundary forcing for the tracers (temperature and salinity) is derived from the hourly Copernicus Marine Environment Monitoring Service (CMEMS) Forecast Ocean Assimilation Model (FOAM) Atlantic Margin Model version 7 (AMM7) output (O'Dea et al., 2012), which is interpolated over the model grid and applied as hourly interpolated vertical profiles at the lateral boundaries using the flow relaxation scheme (FRS) (David, 1976;Engerdahl, 1995). The boundary forcing for water levels and currents is split into three components: a tidal harmonic signal, a barotropic signal and a baroclinic anomaly profile. The tidal harmonic forcing is reconstructed for each model time-step from the tidal constituents for the M2, S2, N2, K2, K1, O1, Q1, P1, and M4 constituents derived from the TPXOv8 model (OSU-OTIS [Egbert & Erofeeva., 2002]). The barotropic forcing consists of the tidally averaged sea surface elevation and depth mean currents and is derived from hourly CMEMS FOAM AMM7 output. The baroclinic forcing is an anomaly of the current profile with respect to the combined tidal and barotropic signal. The Flather radiation scheme (FRS) (Flather, 1994) is used for the tidal harmonic and barotropic forcing at the first lateral boundary bin of the model. The baroclinic forcing uses the FRS. In additional to the lateral tidal forcing, a tidal potential forcing with the same tidal constituents is applied over the whole model domain.
Atmospheric forcing is introduced into the model using the NEMO CORE bulk formulation (Large & Yeager, 2004) and hourly atmospheric forcing fields, which are derived from the European Centre for Medium-Range Weather Forecasts (ECMWF) Reanalysis version 5 (ERA5) data set and consist of 10 m wind components, 2 m temperature and dew-point temperature, mean sea level pressure, and downward solar and thermal radiations (Hersbach et al., 2020). The penetrative solar radiation scheme uses a 2-band approach with the Jerlov water classification type IV parameterization (abs = 0.8, si0 = 0.9, and si1 = 2.10). Surface freshwater input is provided by the ERA5 hourly snowfall and total precipitation in addition to a processed daily river climatology based on river discharge data sets derived from BSHE-HYPE and Met-Oce Figure A1. Sketch of the σ-z* hybrid grid use in NEMO model.

Journal of Geophysical Research: Oceans Appendix C: Estimates of and Issues Regarding Other Balance Operations
In Equation 15, the parameter ρ r and coefficients A, B, and C are given by Millero and Poisso (1981): with these parameters, yields thermal expansion coefficient of seawater at 1 atm for different salinity values ( Figure C1). For the seawater with temperature between 0°C and 30°C, α varies from −0.5 to 3.5 × 10 −4 K −1 . The results are consistent with those calculated by using TEOS-10 (see appendix to Rainer, 2011).
In this study, a 3DVAR system (i.e., GALATON3DVAR) is applied for DA. The most common approach is to use a three-step balance operator (Mogensen et al., 2012), where in the first step, salinity is adjusted based on temperature corrections. Here, water masses (characterized by temperature and salinity) are maintained, and model errors are attributed to vertical positioning errors. This correction only makes sense if there is indeed a vertical salinity gradient in the depth range where temperature is corrected. As the temperature corrections are mainly performed in the mixed water layer, where the salinity gradient is small, such a balance operation is not implemented. The second common balance operation refers to the pressure gradients. In the Boussinesq approximation used for models such as NEMO, the introduction of horizontal temperature gradients by the analysis leads to unrealistic pressure gradients, which can be compensated by a respective balance operator that corrects sea surface elevation. In the current study, such a correction is done by linearly balancing sea surface height with assimilation-changed seawater temperature. The third balancing operation is the adjustment of currents to achieve geostrophic balance of the analysis corrections (Mogensen et al., 2012). However, the assumption of geostrophic balance is unjustified and this approach is not applicable in our case because of the strong influence of bottom friction on the shelf. Figure C1. Thermal expansion coefficient of seawater at 1 atm as a function of temperature and salinities between 0 and 42 psu, as indicated by the curves.