NEMO on the shelf: assessment of the Iberia-Biscay-Ireland configuration

Abstract. This work describes the design and validation of a high-resolution (1/36°) ocean forecasting model over the "Iberian–Biscay–Irish" (IBI) area. The system has been set-up using the NEMO model (Nucleus for European Modelling of the Ocean). New developments have been incorporated in NEMO to make it suitable to open- as well as coastal-ocean modelling. In this paper, we pursue three main objectives: (1) to give an overview of the model configuration used for the simulations; (2) to give a broad-brush account of one particular aspect of this work, namely consistency verification; this type of validation is conducted upstream of the implementation of the system before it is used for production and routinely validated; it is meant to guide model development in identifying gross deficiencies in the modelling of several key physical processes; and (3) to show that such a regional modelling system has potential as a complement to patchy observations (an integrated approach) to give information on non-observed physical quantities and to provide links between observations by identifying broader-scale patterns and processes. We concentrate on the year 2008. We first provide domain-wide consistency verification results in terms of barotropic tides, transports, sea surface temperature and stratification. We then focus on two dynamical subregions: the Celtic shelves and the Bay of Biscay slope and deep regions. The model–data consistency is checked for variables and processes such as tidal currents, tidal fronts, internal tides and residual elevation. We also examine the representation in the model of a seasonal pattern of the Bay of Biscay circulation: the warm extension of the Iberian Poleward Current along the northern Spanish coast (Navidad event) in the winter of 2007–2008.


Introduction
The Northeast Atlantic (NEATL) region covers areas of important economic and social activities that include fisheries, transportation of oil and gas, commercial ship traffic, coastal management, coastal protection and energy production. The availability of validated estimates and forecasts of marine variables in this coastal region is expected to accompany the current development of user-driven activities and applications. The Iberia-Biscay-Ireland (IBI) regional modelling system is one of the seven MyOcean Monitoring and Forecasting Centres (http://www.myocean.eu), and has been developed with those needs in mind. The MyOcean and My-Ocean2 projects aim at developing Marine Core Services within the Copernicus initiative (ex-GMES, Global Monitoring for Environment and Security). The IBI system is to serve downstream coastal modelling users and hence contribute to the development of better products for final users. The operational version of the IBI system described in this paper is now run daily, feeding the MyOcean data access portal with daily analyses and forecasts.
The IBI modelling domain ( Fig. 1) has the singular property of concentrating a large spectrum of physical ocean Published by Copernicus Publications on behalf of the European Geosciences Union. processes, including, for instance, large eddies spawned by the North Atlantic Current and the Azores Current, intense upwellings along Portuguese coasts and gravity currents flowing down the Gibraltar and Faroe channels. Characteristic to the region is a relatively steep slope separating the deep ocean from the shelf. Along these, poleward slope currents flow in the subsurface; they are observed as far north as Ireland (White and Bowyer, 1997). In the presence of stratification, strong bathymetric gradients trigger in very localized spots the conversion of barotropic tides into internal waves which propagate onto the shelf and into the deep ocean (Pingree et al., 1986). On the continental shelves, intense tidal motions provide the dominant source of energy. The associated turbulent-mixing shapes much of the water column properties, preventing under some conditions any surface stratification to setup (Simpson and Hunter, 1974).
From a modelling point of view, this remarkable variety of processes and scales is, however, particularly challenging. Historically, coastal and deep-ocean numerical models have followed more or less separate paths. The NEMO model (Nucleus for European Modelling of the Ocean; Madec, 2008) used in this study has been developed essentially to perform global climatic simulations. The recent desire shared by several European operational oceanography centres to hold a single tool for both global and regional applications has strongly influenced its development to that end (O'Dea et al., 2012). It is likely that, in turn, global-scale simulations will benefit from more realistic coastal modelling performance. As already explored in other models, explicit representation of tidal motions is likely to become a standard feature of global eddying models in the near future (see Arbic et al., 2010 for instance). Also, one particular aspect in the present model configuration is the relatively high (2-3 km) horizontal resolution used. This clearly pushes the model into the sub-mesoscale-permitting regime over much of the domain and allows a significant part of the internal wave spectrum to be resolved.
In this paper, we attempt to give a broad-brush account of one particular aspect of setting up an ocean model, namely consistency verification. This type of validation is often conducted upstream of the implementation of the system which will later be used for production and routinely validated; it is conducted by comparing model results with observations, and is meant to guide model development in identifying gross deficiencies in the modelling of several key physical processes. This work is applied research that develops a scientific framework and methodology for improving ocean model configurations at the development stage for use in basic research or operations. This paper's approach is partly inspired by the works of , Holt et al. ( , 2005 and Sotillo et al. (2007), partly by the specific needs of this project and partly by available observational data for this project. In particular, a challenging issue is the consistency verification of high-frequency dynamics, such as barotropic and internal tides, or circulation variability at daily timescales. Additionally, we aim at illustrating the potential use of this regional system to provide ocean state estimates which can help interpolate between local or patchy observations and integrate them in larger-scale patterns and processes. To better collaborate and meet My Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ Ocean objectives, we focus on evaluating the model configuration during the year 2008. The paper is organized as follows. Section 2 focuses on the NEMO model, its new developments and the experimental protocol. Section 3 provides domain-wide consistency verification results in terms of transports, barotropic tides, sea surface temperature (SST) and upper layer stratification. Section 4 focuses on two dynamical subregions: the Celtic shelves and the Bay of Biscay slope and deep regions. The model-data consistency is checked for variables and processes such as tidal currents, tidal fronts, internal tides and residual elevation. We also examine the winter warm extension of the Iberian Poleward Current along the northern Spanish coast (Navidad event) in the winter of 2007-2008. Finally, Sect. 5 gives a summary, a discussion and perspectives for improving the model.

NEMO physics and numerical aspects
The IBI model numerical core is based on the NEMO v2.3 ocean general circulation model (Madec et al., 1998;Madec, 2008). It solves the three-dimensional primitive equations in spherical coordinates, discretized on an Arakawa C-grid, assuming hydrostatic and Boussinesq approximations. NEMO is written in a generalized vertical coordinate framework allowing for some flexibility in the choice of the vertical discretization. Hence, one can easily switch from purely geopotential coordinates (with optionally partial bottom cells as used here) to, for instance, terrain-following vertical coordinates. The bulk of the numerical code used here is very similar to the one described in Barnier et al. (2006). This includes the vector invariant form of the momentum equations and the energy-enstrophy discretization of vorticity terms. There are, however, important differences related to the specific coastal and tidal dynamics studied here that are briefly described below.
In order to properly simulate tidal waves, the default "filtered" free-surface formulation (Roullet and Madec, 2000) must be discarded. This scheme has indeed been designed to damp fast external gravity waves in order to extend the permissible time step. As an alternative, a conventional timesplitting scheme is used here: the barotropic part of the dynamical equations is integrated explicitly with a short time step (3 s), while the more costly update of depth-varying prognostic variables (baroclinic velocities and tracers) is carried out with a larger time step (150 s). The mode-coupling procedure as well as the barotropic time-stepping scheme (a "generalized forward backward") follows the work of Shchepetkin and McWilliams (2004) adapted to the standard NEMO leap-frog time stepping of depth-varying variables.
Because of the explicit simulations of tides, sea level elevation can become large on the shelf compared to the local depth with the result that linear free-surface approximation needs to be relaxed. In practice, all model vertical thicknesses dz are remapped in the vertical at each time step to account for the varying fluid height according to dz(x, y, z, t) = dz 0 (z) H (x, y) + η(x, y, t) H (x, y) . (1) Here dz 0 stands for the "reference" vertical thicknesses, η is the sea level, H is the water depth at rest, (x, y, z) the model geographical coordinates and t the model time step. This vertical coordinate system is called z * coordinate (Adcroft and Campin, 2004). This implies a correction (here in a standard density Jacobian form) of the hydrostatic pressure gradient force since computational surfaces are not horizontal anymore. The use of z * coordinates reduces computational error significantly compared to terrain-following coordinates, with z * surfaces having very faint slopes (Marsaleix et al., 2009). The model turbulent-mixing scheme uses parameterization and equations from Warner et al. (2005) unless mentioned explicitly here. Vertical turbulent-mixing processes are parameterized with a k-epsilon two-equation model implemented in the generic form proposed by Umlauf and Burchard (2002). The model is complemented with the type "A" full equilibrium form of Canuto et al. (2001) stability functions. The turbulent kinetic energy background that accounts for unresolved internal wave breaking is set to k min = 10 −6 m 2 s −2 . Dissipation is limited under stable stratification using a Galerpin coefficient c galp = 0.267. As discussed by Holt and Umlauf (2008), these two important parameters control the magnitude of the background viscosities/diffusivities away from boundary layers.
Flux boundary conditions are used for their better accuracy (Burchard et al., 2005). Compared to clamped boundary conditions, Burchard et al. (2005) show that they better handle changes in vertical resolution, which is an important requirement here with the highly variable partial bottom cell thicknesses ( Fig. 2a). At the bottom a steady balance between shear and dissipation characteristic of log-law behaviour is assumed. At the surface, turbulent kinetic injection through surface wave breaking is considered leading to a surface kinetic energy flux (Craig and Banner, 1994): where α =100 is a parameter and u * the ocean surface friction velocity scale deduced from the wind stress. Related to the wave breaking parameterization is the still largely uncertain choice for the surface roughness zo s that brings the final piece in the vertical-mixing boundary conditions. A common assumption is to scale it with the wind sea wave height H s : γ is a free parameter ranging from 0.6 to 1.6 in the literature (γ = 1.3 in our case). If no observational nor modelling www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 data are available, various empirical formulae express H s as a function of wind stress. We followed here the formulation of Rascle et al. (2008) so that where g is the gravity acceleration, w * the atmospheric friction velocity and w * ref =0.6 m s −1 a typical wind velocity scale above which wave growth is limited. Rascle et al. (2008) showed that Eq. (4b) better fits their model dataset compared to the classical approach that sets β to a constant value. From the various experiments we made, use of this formulation effectively tapers the effect of wave breaking in the case of strong winds, which was found to substantially improve the sea surface temperature (comparison with surface buoys). Along lateral boundaries, free-slip boundary conditions are used everywhere except inside the Strait of Gibraltar, where no slip conditions are applied. This prevents the Alboran jet in the Mediterranean Sea from turning into an unrealistic quasi-permanent "coastal mode" (i.e. trapped along the African coast). At the bottom, a quadratic bottom drag with a logarithmic formulation is calculated as where κ = 0.4 is the Von Karman constant, Cd min = 2.5 × 10 −3 a minimum drag coefficient, dz b is the lowermost bottom cell thickness and z 0b = 3.5 × 10 −3 m is the bottom roughness. With the reference geopotential discretization chosen here, the logarithmic formulation applies for depths shallower than 170 m, i.e. on most of the continental shelf (Fig. 2b). Although ocean models interchangeably assume constant or logarithmic formulations as bottom drag formulation, the use of Eq. (5) is mandatory here to ensure a continuous (i.e. not related to numerical settings such as the vertical discretization) representation of bottom stress as it should be in nature. This somewhat counterintuitive property owing to the largely spatially discontinuous drag coefficient ( Fig. 2b) comes from the tight interplay between the modelled vertical shear and the vertical-mixing scheme (impact on model stability of this interplay was further evidenced by Burchard et al., 2005). The formulation of the bottom stress must agree with the assumed law of the wall in the verticalmixing scheme, a requirement emphasized here by the use of partial bottom cells.
Tracers' advection is computed with the QUICKEST scheme developed by Leonard (1979). This third-order scheme is well suited to the high resolution used here and the modelling of the sharp front characteristics of coastal environments. Lateral sub-grid-scale mixing are parameterized according to horizontal biharmonic operators for both momentum (A m = −2.5×10 8 m 4 s −1 ) and tracers (A t = −2.5× 10 7 m 4 s −1 ). Note that the latter value is particularly small, since significant inherent diffusion is also associated with the QUICKEST scheme.

Model setup
The model domain covers the Northeast Atlantic Ocean from the Canary Islands to Iceland and from about 20 • W to 10 • E (Fig. 1). Even if the focus of the IBI forecasting system remains on the Iberian, Biscay and Irish regions, the computational domain has been extended to include the western Mediterranean Sea and the North Sea. This strategy has been chosen to lessen the impact of northern and eastern boundaries with the parent model, PSY2V3. Indeed the PSY2V3 model, an operational forecasting system covering the North Atlantic and the Mediterranean Sea with a 1/12 • resolution, does not simulate all the physical processes of interest here (tides in particular). Moreover, such a strategy allows for the explicit computation of mass, heat and salt fluxes in the narrow Gibraltar and Kattegat straits, which we expect will be better resolved thanks to the improved horizontal resolution.
The primitive equations are discretized on an horizontal curvilinear grid which is a refined subset at 1/36 • (≈ 2-3 km) of the so-called "ORCA" tripolar grid, commonly used in other NEMO-based large-scale and global modelling experiments (Barnier et al., 2006). This is also an exact 3 : 1 refinement of the PSY2V3 model grid that provides initial and boundary conditions. This strategy greatly simplifies interpolation procedures and allows for exact boundary fluxes/bathymetry matching, the latter increasing the overall robustness of the coupling. If we consider the common rule of thumb that sets the internal Rossby radius as the minimum horizontal grid spacing to resolve the first baroclinic instability mode, the chosen grid clearly moves the model into the eddy-resolving regime, at least over the deep ocean (the minimum Rossby radius as computed from climatology and for depths greater than 1000 m is around 5 km). On the shelf, much lower internal radii (2 km) only allow the model to be in an "eddy-permitting" regime. The same 50 reference z levels as in the parent grid model are used in the vertical, with a resolution decreasing from ∼1 m in the upper 10 m to more than 400 m in the deep ocean. A partial step representation of the very last bottom wet cell is used with some constraints on the resulting minimum bottom cell thickness to guarantee model stability (it must be greater than 15 m or 20 % of the reference grid thickness).
The original bathymetry is derived from the 30 arc-second resolution GEBCO 08 dataset (Becker et al., 2009) merged with several local databases (F. Lyard, personal communication, 2010). The model bathymetry is linearly interpolated from the resulting composite and slightly smoothed using a Shapiro filter to remove grid scale noise and remaining discontinuities between local databases. At open boundaries, Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ within 30-point-wide relaxation areas, bathymetry is exactly set to the parent grid model bathymetry and progressively merged with the interpolated dataset described above. This approach and the identical vertical grid as in the parent system imply that no vertical extrapolation of open boundary data is necessary. Since no wetting and drying capability is presently available in NEMO, a minimum value for H , H min , has to be chosen to ensure that the fluid height h = H + η always remains positive in Eq. (1). Sensitivity tests revealed that tidal propagation in the North Sea was highly dependent on the chosen H min value if simply set, for instance, to a uniform value greater than the maximum tidal amplitude in the area (≈ 10 m). We finally set where Atide max (x, y) is the maximum tidal amplitude (we use here as a rough estimate the sum of FES2004 elevation amplitude harmonics; see Lyard et al., 2006, for a description of the FES2004 tidal atlas).

Atmospheric forcing
Meteorological fields from the European Centre for Medium-Range Weather Forecasts (ECMWF) with a 3 h period and 0.25 • horizontal resolution drive the present model simulations. According to Bernie et al. (2005), this temporal resolution associated with adequate vertical resolution in surface layers (typically 1 m, as used here, or less) is sufficient to model diurnal variations of SST. Evaporation, latent and sensible heat fluxes, and wind stresses are computed according to Large and Yeager (2004) bulk formulae. Net longwave heat flux is obtained from ECMWF downward flux to which the following blackbody radiation is added: where σ is the Stefan-Boltzman constant and E = 0.99 the emissivity.
Penetrative solar radiation Q sr is parameterized according to a two-band exponential scheme such that www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 The first exponential term from the left corresponds to red and near-infrared radiation absorbed in surface layers (the chosen e-folding length scale is l = 0.35 m), while the second refers to shorter wavelengths mostly in visible and ultraviolet bands. R = 0.54 determines the fraction in each band of the available penetrative solar radiation at the surface Q sr (0) (we assume a constant albedo of 6.6 %). The latter part of the spectrum penetrates at much greater depths and participates to photosynthesis such that it is often referred to as photosynthetically active radiation (PAR). While there is little variation of absorption depth in the first band, PAR absorption coefficient (k PAR ) is highly dependent on ocean turbidity. Since the IBI model domain encompasses areas with very different optical water properties (from highly turbid on the shelf to clear waters in the subtropical gyre; see Fig. 2c), a spatially variable k PAR coefficient has been used. It is built from a 10 yr, 9 km, monthly climatology of SeaWiFS satellite diffusive attenuation coefficient at 490 nm (kd 490 ), transformed into an equivalent PAR attenuation coefficient according to Morel et al. (2007): Equation (9) and the different algorithms used for the retrieval of attenuation from satellite data are nevertheless only valid for Case-1 waters (where chlorophyll concentration controls optical properties). As a correction, Seawifs-based estimation is merged with a monthly climatology processed by Ifremer (Gohin et al., 2005) which is valid in coastal waters. The k PAR threshold value between Case-1 and Case-2 waters (over which only Ifremer data are considered) has been set to 0.12 m −1 . This value is more or less the lower bound of agreement of the Gohin et al. (2005) estimate if compared to in situ measurements.

River runoffs
Climatological monthly flow rates are prescribed for 33 rivermouth locations (Fig. 2c). These have been obtained by averaging data from the Global Runoff Data Centre (http://grdc. bafg.de) and the French hydrographic database "Banque Hydro" (http://hydro.eaufrance.fr). Rivers are applied by specifying (i) a constant velocity in the vertical whose integral matches the specified transport, (ii) Neumann condition for temperature and (iii) a constant salinity (0.1 psu).

Open boundaries
After numerous tests, relatively simple yet robust open boundary formulations have been chosen. The method of characteristics is used for barotropic variables following Blayo and Debreu (2005). These conditions are expressed (for a western open boundary) as where U and V refer respectively to normal and tangential barotropic velocities; g is the gravity acceleration, h the water depth and η the sea level. Index "b" corresponds to the prescribed boundary value, index "i" to the model-computed value immediately inside the domain and index "ext" to external data. A Neumann condition for elevation is used, while baroclinic velocities, temperature and salinity are prescribed to external data. For the latter variables a 30-point relaxation area with a minimum timescale of 1 day is used, which strongly damps outgoing perturbations from the prescribed fields, while allowing for the explicit simulation of highfrequency fluctuations related to the atmosphere in the surface layers. Temperature, salinity, velocities and sea surface height (hereafter SSH) from PSY2V3 daily outputs are used as the slow component of open boundary data. PSY2V3 does include a multivariate, Kalman-based weekly data assimilation procedure of along-track altimetry data, SST and hydrological observations (Drévillon et al., 2008;Dombrowsky et al., 2009). SSH and barotropic velocities tidal components are then added as the sum of a maximum of 11 constituents (M2, S2, K2, N2, K1, O1, P1, Q1, M4, Mf, Mm) provided by the TPXO 7.1 global tide model (Egbert et al., 1994). A 35-day barotropic experiment (i.e. without stratification), with clamped boundary conditions for sea level and Neumann conditions for velocities, is computed to estimate the 6 major barotropic tidal velocities constituents used for the open boundary conditions (similar to the runbt experiment described below). This method, reminiscent of the iterative procedure of Flather (1987), clearly improved the overall tidal statistics (semidiurnal phases otherwise exhibit significant errors due to the inadequacy of TPXO velocities with the model physics).
Finally, since surface atmospheric pressure forcing is explicitly considered in the dynamical equations, the inverse barometer signal (IB) is added to sea level open boundary data. In the Atlantic part of the domain, IB is a fairly good approximation of the sea level response to pressure variations, while it has serious limitations in a semi-enclosed sea like the Mediterranean Sea (Le Traon and Gauzelin, 1997). Since we essentially focus on the Atlantic side of the basin, we have left for further study the use of more-sophisticated methods (such as the Candela (1991) analytical model in the Mediterranean Sea) or outputs from Baltic and Mediterranean Sea models that do explicitly include the pressure component. As demonstrated in the following, this simplification still allows for pretty realistic variations of Mediterranean outflow transport at Gibraltar, which can be considered as the boundary of our domain of interest.

Experiments
Results from two simulations will be considered in the next sections: 1. A "tide-only" run (hereafter "runbt") which consists in a 35-day "barotropic-like" simulation (homogenous density, no atmospheric forcing, only tidal forcing). Since proper representation of bottom boundary layers is also of interest here, this experiment retains the 50 vertical levels distribution and the k-epsilon vertical-mixing scheme described above. Compared to realistic simulations, it will provide insights on the impact of stratification, and, in particular, changes due to internal waves on the simulation of tidal waves. Note that a reduced set of seven harmonics (M2, S2, N2, K1, O1, Q1, M4) is used since these waves are separable through harmonic analysis over such a short simulation period. Since these harmonics represent most of the tidal energy, it is expected that this will not hamper the comparison to the realistic experiment.
2. A realistic simulation with full forcing and stratification (hereafter "runbc"). It has been initialized on July 25 2007 from temperature, salinity, sea level and velocities linearly interpolated from PSY2V3 analysis. Tidal and atmospheric pressure forcing are smoothly introduced (both at open boundaries and in dynamical equations) thanks to a 2-day linear ramp. Meanwhile, horizontal viscosity linearly decreases from 4 times its nominal value to damp instabilities arising from the use of unbalanced initial fields. The simulation ends in February 2009. In the next sections, we focus on the year 2008 only.

General comments
Since tidal motions are the dominant source of energy on the continental shelf, we give here some quantitative assessment of model performance on that particular aspect. Numerous modelling studies have already conducted such an exercise (Davies and Aldridge, 1993;Davies et al., 1997;, which has proven to be sensitive to various numerical details, due in particular to the complex tidal propagation pattern in the English Channel and the North Sea.
In order to keep the present account synthetic enough, a detailed description of the different tests performed to reach the present results cannot be given. The use of "equilibrated" boundary conditions as described in Sect. 2.3.3 is obviously one of the keys as well as the use of the tide-generating force owing to the size of the domain.
Here we solve this problem with the use of geopotential coordinates. Terrain-following coordinates have historically been preferred to maintain sufficient vertical resolution at the bottom and to have a continuous representation of the bathymetry. The use of partial bottom cells nevertheless partially restores the latter property with z coordinates. Moreover, the vertical discretization chosen here sets the level of the last wet cell at most at 10 m above the bottom for depths lower than 150 m (Fig. 2a), which should be enough to properly resolve the bottom boundary layer.
An aspect that also deserves some attention is the impact of stratification on tidal amplitudes. The high resolution used here is expected to allow for a significant part of the tidally induced internal wave spectrum to be resolved, which in turn provides an important sink of energy for barotropic tides. In forward global tide modelling, correctly accounting for this barotropic to baroclinic conversion explicitly or through adhoc parameterizations has now become key for further improvement (Arbic et al., 2010;Lyard et al., 2006). In the next paragraphs, we examine results of a 1 yr harmonic analysis on the runbc experiment. These are eventually compared to the barotropic experiment runbt with harmonic analysis performed over the last 30 days.

Tidal elevations
The M2 co-tidal chart is shown in Fig. 3a and depicts the complex wave propagation around amphidromes in the North Sea. The overall picture is in good agreement with published charts, including the position of the degenerate amphidrome south of Norway that has proven to be difficult to obtain in other modelling studies (O'Dea et al., 2012). Further insights of the model performance can be gained by comparing the result to FES2004 atlases (Lyard et al., 2006). The FES2004 solution has been obtained through the use of a finite-element hydrodynamical model assimilating altimetry and in situ harmonic data; it can be considered a fairly reliable reference, at least in the deep ocean. Figure 3b shows M2 difference with the FES2004 solution, amplitude differences being generally smaller than 1 cm in the deep ocean (z > 1000 m), phase errors (not shown) being negligible. Most of the differences there are related to the surface signature of internal waves, which are further investigated in Sect. 4.3.2. On the shelf, the regional model is close to the FES2004 solution, except in the Dover Strait area and in the German Bight. From the same comparison performed in the barotropic experiment, it is interesting to note that much larger differences arise when no stratification effects are present (Fig. 3c). The M2 wave amplitude in that case is systematically overestimated, a difference reaching more than 15 cm at the Western English Channel entrance and in the Irish Sea (phase is not significantly different between experiments). From the pattern of the amplitude difference between experiments (not shown), it is tempting to relate this difference to barotropic to baroclinic energy transfer: wave amplitude indeed undergoes a www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 sharp decrease at 48 • N along a well-observed spot of internal wave generation. Since other factors may explain the differences between experiments (effect of stratification on viscosity profiles for instance), a separate study beyond the present "global" assessment would, however, be needed to firmly reach that conclusion.
To gain some more quantitative measure of errors, comparisons are made at various tide gauge locations for the four main tidal constituents (M 2 , S 2 , K 1 , O 1 ). Comparison for M 4 compound wave is also considered, as an indicator of how well non-linear interactions are modelled, in particular those induced by the z * coordinate (Eq. 1). For each tidal component, the RMS of the complex difference (see Davies et al., 1997 for the definition), which synthesizes both amplitude and phase errors, is given in Table 1. To appreciate how significant these errors are, it is instructive to make a comparison with modelling studies based on well-validated shelf sea models such as the POLCOMS model used in Holt et al. ( , 2005. Considering the overlapping area between both models, the present experiment has similar accuracy for M2 (21.6 cm compared to 24.1 cm in IBI). IBI has larger errors for the M 4 component, but it seems to outperform the  model on the other components, with errors reduced by more than a half for the diurnal component. We however stress that measurement points differ between our study and that of  the present dataset being essentially distributed along the coast.

Barotropic tidal currents
Only a brief account is given here on barotropic tidal velocities since a more detailed discussion on total surface tidal currents is given in the regional assessment section. Histori-cal current moorings collected within the framework of the WOCE experiment are used for comparisons of modelled barotropic tidal ellipses (Fig. 1, orange dots). Table 2 gathers comparisons to in situ observations shown in Fig. 4 for the dominant M2 tidal component in the runbt simulation. The overall agreement is rather good with a global RMS error of 4.5 cm s −1 for the semi-major axis, model velocities being generally slightly larger than observations. Along the shelf slope in the northern part of the Bay of Biscay, where plenty of measurement are available, absolute errors are smaller than 8 cm s −1 (75 % are smaller than 4 cm s −1 ). This is of similar accuracy to the regional model of Pairaud et al. (2010). Results with full forcing and stratification exhibit very similar errors (RMS = 3.8 cm s −1 ), velocities being, however, slightly underestimated.

Transports
Transport estimates give access to broad lines of the model circulation. Monthly volume, heat and freshwater transports are computed across various sections and are compared to previous published estimates ( Table 3). Note that published transport values from observations do not always contain transport error estimates. The standard deviation of published estimates from the monthly model transports is used to give an order of magnitude of the errors between model and observations (Table 3). This value includes model errors as well as transport interannual variability.
For Sect. 1, the transports associated with the Azores Current and the Mediterranean Water (MW) outflow (15.1 Sv and −12.3 Sv) are larger than the eastward and westward Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ observed transports of 13.7 Sv and −4.7 Sv (Peliz et al., 2007). However, the observed estimates are deduced from data obtained in September-October of two consecutive years in 1991 and 1992, and the modelled estimates show high monthly variability, indicating that the eastward flow is of comparable magnitude with observations, while the westward flow is overestimated. Across the Strait of Gibraltar, model transports have been computed from monthly averaged values of currents and salinity (Atlantic Water inflow of 0.48 Sv and the MW outflow of −0.49 Sv). Tsimplis and Bryden (2000) found values of 0.46 Sv and −0.35 Sv at the Strait of Gibraltar using currents averaged over several months (estimations not included in Table 3). However, the flow variability is dominated by high-frequency tidal currents through the Strait of Gibraltar (Tsimplis and Bryden, 2000) and transports across the Strait of Gibraltar are commonly computed from hourly or daily fields rather than from yearly averaged salinity and velocity fields. Using this method, we obtain transports of 0.51 Sv in the upper layer and of −0.71 Sv in the lower layer, which lie within the range of observation-based estimates. We also note that the outflow is larger than the inflow, while the inflow is normally larger to compensate for the excess of evaporation minus precipitation over the Mediterranean Sea. In the IBI solution, the outflow transport is mainly influenced by the open boundary conditions prescribed from the PSY2V3 solution which explains the larger outflow in the strait. Table 2. RMS between observed and modelled (runbt) tidal ellipse parameters. Units are in cm s −1 for the semi-major axis (SEMA) and the semi-minor axis (SEMI) and in • for the inclination (Inc) and the phase (Pha). Statistics are also performed over the area covered by the studies of Holt et al. ( , 2005 , and results obtained in the studies of Holt et al. ( , 2005 Mauritzen et al., 2001). Previous estimates use a density of 31.7 to separate inflow of Atlantic Water from outflow of MW, but comparisons of hydrological fields with climatology have shown than the MW is too light in the model (not shown) and the 31.7 isodensity may not be appropriate to separate the two water masses. When choosing a density threshold value of 30.3 in the Gulf of Cadiz region, we obtain transports of 2.0 Sv and −1.1 Sv across the Gulf of Cadiz and a transport of −1.3 Sv for the MW. The results are closer to the estimates from observations and the transport lies within the range of published values for the inflow, but the westward flow is still underestimated (Howe, 1984;Rhein and Hinrichsen, 1993;Mauritzen et al., 2001). South of Portugal, the MW flows northward with a transport comparable with previous estimations (Da Silva, 1996;Mazé et al., 1997;Coelho et al., 2002). Off the west of the Iberian Peninsula the transport is slightly overestimated at 40 • N but is in good agreement with previous estimates at 43 1996; Mazé et al., 1997;Mauritzen et al., 2001;Coelho et al., 2002;Alvarez et al., 2004;Lherminier et al., 2007).
Transports of the slope current across the Celtic shelf slope and along the Ellett line are of the same orders of magnitude than published values (Ellett and Martin, 1973;Pingree and Le Cann, 1990;Holliday et al., 2000). However, we note that the transport across the Celtic slope is mainly constrained within the upper slope in the model, while it is more evenly distributed along the slope in the study of Pingree and Le Cann (1990). The distribution of the transport with depth is comparable with estimations from data across the Ellett line section. Through the Dover Strait, both observed and modelled transports are weak (Otto et al., 1990;Prandle et al., 1996), but the transport is slightly underestimated in IBI.

Mediterranean outflow
The transport through the Strait of Gibraltar is controlled by various processes at different time scales. Among these processes, tides are the most energetic one influencing the mean flow (Tsimplis and Bryden, 2000). Thus hourly trans-port through the Strait of Gibraltar has been calculated to examine the outflow variability. Following the study of Sánchez Román et al. (2009) the interface to separate the inflow from the outflow has been defined as the time-dependent depth of the surface of zero low-pass frequency velocity. The modelled outflow transport estimated with this interface is denoted TMV (transport model velocity). The outflow transport computed using the commonly used 37.25 isohaline (TMS, transport model salinity) is also computed for comparisons. Figure 5 illustrates time series of observed and modelled transports for TMV and TMS over a period overlapping the model simulation and the transport estimates of Sánchez Román et al. (2009), denoted TOV (observation velocity). The definition of the interface based on the depth of the zero velocity gives better results than the isohaline interface. In comparison to the TOV mean transport of −0.72 Sv, the mean outflow of TMV is −0.74 Sv and the correlation coefficient with TOV is 0.71, while a lesser agreement is found for TMS with a too small transport of −0.69 Sv and a correlation coefficient of 0.66. Differences between peak values and lower values at neap tides reach 0.70 Sv both in TOV and TMV.

SST datasets
We use SST processed at Météo-France/CMS for comparisons with the model SST (Le Borgne et al., 2011). Data consist in L3 multi-sensor products built from bias-corrected L3 mono-sensor products. SST satellite data fields are produced daily with a 0.05 • resolution and an accuracy of 0.5 • C. Sea surface temperature measured from in situ buoys is also used (Fig. 1, Table 4). The dataset consists in hourly time series which enable the computation of the diurnal cycle. Figure 6 represents maps of seasonal SST biases and RMS of differences between the satellite data and the model fields, for summer (July-August-September) and winter (January-February-March). The mean bias over the whole IBI area is low (less than 0.25 • C) over the whole year but is not evenly distributed in space and time. In winter, the bias is lower than 0.2 • C and the RMS of the difference is on average below 0.5 • C. Discrepancies are higher in summer, with biases locally larger than 1 • C over the Armorican shelf, in the English Channel and west of Iberia. West of Iberia the model is too warm in both seasons, with the SST overestimation being larger in summer during upwelling season. In the Irish Sea, at the entrance of the English Channel and in the Iroise Sea, the discrepancies of 1-1.5 • C may be explained by the high variability of some thermal front positions which are tidally induced, as discussed in Sect. 4.2. In the Gulf of Cadiz (GC) the modelled SST is too warm along the slope, which reveals Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/  20N Prandle et al. (1996) the model's relative failure at representing the GC slope current (GCC). The GCC, intimately linked to the westward MW outflow at depth (Peliz et al., 2009), indeed advects cold upwelled water along the western Portuguese coast into the GC. This stresses the model's limits, already outlined above, at representing MW outflow in its present version. Figure 7a represents a Taylor diagram of the observed and modelled SST at buoy locations (Fig. 1, Table 4). SST time series have been previously 1-day low-pass-filtered to remove the diurnal cycle signal, which is discussed in the next section. SST is particularly well modelled in IBI with "normalized" standard deviation close to 1 and correlation greater than 0.95 (except at Cabo Silleiro and Villano-Sisargas, www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 where the correlation coefficients are 0.85 and 0.89, respectively; these two locations are included in the "Iberia" subregion described in Fig. 7). When looking regionally, larger discrepancies are found along the Armorican slope (included in the "Biscay" region) and within the Gulf of Cadiz (included in the Iberia region) due to overestimated modelled SST during summer. Warmer model SST during summer at these locations leads to a larger seasonal amplitude. The IBI model ability to reproduce the SST variability results from the high spatial resolution and from model developments (see Sect. 2); the latter allow the model to simulate and predict slope and shelf processes (e.g. tidal mixing) better than deep-ocean configurations such as PSY2V3 (however, SST data assimilation leads to closer agreement with observations in PSY2V3). The model-data discrepancies in the upwelling areas suggest a need for investigation of the sensitivity of the model response in these areas to the wind product used to force the model.

Diurnal SST cycle
The combination of daily variations of atmospheric fluxes and weak winds can lead to strong SST diurnal anomalies (hereafter Dsst). The Dsst amplitudes are computed over a 24 h period as the difference between the maximum SST and the minimum SST previously corrected from the daily trend (Pimentel et al., 2008). The Dsst for the buoy data are only computed when the 24 SST values per day are available to avoid any miscalculations of the Dsst amplitudes. Figure 7b displays box plots of the observed and the modelled Dsst amplitudes over the area covering the Bay of Biscay, the Celtic shelves and the North Sea (areas shown in Fig. 3); it illustrates the distribution of data by showing different values of Dsst amplitudes below which a certain percentage of data are observed. Globally Dsst amplitudes are underestimated in the model, and the median observed Dsst is 0.25 • C, while the median modelled Dsst is only 0.20 • C. When looking regionally at the distribution (not shown), we found that Dsst is underestimated almost everywhere for the 50 % of the data with smaller amplitudes (i.e. Dsst amplitudes below the median value). The resolution of the model and of the forcing, the surface heat fluxes and the vertical mixing in the upper ocean are critical aspects for the diurnal cycle modelling (Pimentel et al., 2008). According to Bernie et al. (2005), a minimum of a 1 m vertical resolution in the upper ocean and of a 3 h temporal resolution of surface fluxes is required to simulate 90 % of the observed Dsst. IBI is forced with 3 h atmospheric fields, but has a coarser Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ vertical resolution with 8 levels (18 levels) in the top 10 m (50 m), which only enables representation of 70 % of the Dsst (Bernie et al., 2005). This partly explains the underestimation of the Dsst in the model. Overestimated wind stress may also lead to an underestimation of the Dsst by inducing too strong vertical turbulent mixing and preventing diurnal heating. As for lower frequencies (see previous section), we recommend further assessment on the realism of the wind field forcing.

Upper layer stratification
The upper layers stratification is a key process to represent in a regional model for its implications on circulation and energetics, in particular at the air-sea interface, and for its strong influence, or even coupling, with biological processes. The subsurface ocean is sparsely sampled with observations. In the Bay of Biscay, for instance, only two moorings of hourly temperature and salinity (T /S) profiles are available to the community (see Sect. 4.3.3). The UK MetOffice provides the so-called EN3 dataset that gathers subsurface T /S profiles from different sensors (e.g. ARGO floats) with quality control based on a comprehensive set of checks (Ingleby and Huddleston, 2007). Figure 8 illustrates the available data for February and August 2008. The sampling is inadequate for a detailed investigation of the processes underlying the T /S distribution. Instead it allows an overall view of the model performance in different areas. Figure 8 shows the RMS of the model-data misfits in temperature over the upper 200 m. We also compute the mixed layer depth (MLD) from the model and EN3 profiles using a criterion on density and a threshold of 0.02 kg m −3 (not shown). In February, modelled temperatures show very good agreement with data almost everywhere, with RMS globally lower than 0.5 • C. The profile of temperature RMS computed over the whole IBI domain (not shown) indicates that maximum RMS values are reached between 50 and 130 m, which is consistent with differences in modelled and observed MLD at these depth ranges. In August, larger RMS in temperature profiles (up to 2 • C) are found on shelves as well as in the deep ocean, with a model overestimation of the temperature in the upper 200 m, in the Iroise Sea and in the English Channel. The MLD is underestimated in the model in summer.
Besides the difficulty of single-point comparison (in particular because of the mesoscale signals), these comparisons emphasize the need to further test the vertical-mixing physics: parameterization and values of parameters, as well as forcing fields (impact of over-or underestimation of wind forcing on turbulence). Another possible source of error may be the wave-mixing parameterization coefficient. Sensitivity experiments performed with a 1/12 • resolution configuration have shown that the modelled MLD (and consequently temperatures in the upper ocean layers) is very sensitive to the wave-mixing parameterization coefficient. The model MLD representation can be enhanced by tuning this parameter, but the parameterization seems to be as inadequate in semienclosed seas as in the Mediterranean Sea or the North Sea. Thus the parameter used is a compromise to minimize the MLD errors both in the open ocean and in the semi-enclosed seas.

Introduction
We now focus on the Bay of Biscay and explore the modeldata consistency for processes specific to the shelf (Sect. 4.2) as well as slope and deep regions (Sect. 4.3). The surface circulation in the Bay of Biscay is relatively weak over the abyssal plain (a few cm s −1 , Charria et al., 2011); it is mainly cyclonic in winter and anticyclonic in summer. Slope currents have been evidenced off the northern Iberian coast and over the Armorican shelf break (Pingree and Le Cann, 1990) Iberian Poleward Current). The circulation over the shelves is dominated by tidal currents; the residual currents are mostly wind-driven and are highly variable in time. Anticyclonic and cyclonic eddies are observed in the deep region; they mostly result from instabilities of the slope current developing at capes and canyons (Pingree and Le Cann, 1992). Internal tides are generated in regions of interaction between barotropic tidal currents and the shelf break. Two sites of generation have been described in the literature: the Armorican shelf break around 47 • N (Pairaud et al., 2010) and the Galician coast (Pichon and Corréard, 2006). We investigate the representation of these processes based upon the available observations. The emphasis is on shelves and slopes because of potential implication for users in biology or other downstream applications of the IBI system. We also examine the representation of the Iberian Poleward Current flow along the northern Spanish coasts in winter. The objective is to show that such a regional modelling system has potential as a complement to observations (an integrated approach) to give information on non-observed physical quantities and to provide links between observations.

Surface tidal currents in the Iroise Sea
In an attempt to locally investigate the performance of the model at reproducing tidal surface ellipses in one particular region characterized by intense tidal currents (up to 4 m s −1 in some interisland passages; SHOM, 1994), we use surface current data from the HF WERA (Wellen Radar) radars at Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ the Brezellec headland, 48 • 04 N,4 • 40 W, and at the Garchine headland, 48 • 30 N,4 • 46 W (Fig. 1). This system is operated by SHOM (the French Navy Hydrographic and Oceanographic Service) and provides radial surface currents every 20 min with a resolution of 1.5 km and an accuracy of a few centimetres per second (Le Boyer et al., 2009). The resulting radial currents have been interpolated on a regular grid with a 2 km spatial resolution extending from 6.78 • W to 4.65 • W and from 47.30 • N to 49.26 • N (Muller et al., 2009) using open-boundary modal analysis (OMA; Kaplan and Lekien, 2007). According to Le Boyer et al. (2009), radar-derived currents have an upper-bound uncertainty of 15 cm s −1 , which still allows us to compute the tidal current components as the latter are larger than the uncertainty and dominate the signal in the Iroise Sea. We use a 6-month dataset in 2007 to perform the harmonic analysis of nine tidal components (M2, S2, K2, N2, K1, O1, P1, Q1, M4). Figure 9a (Fig. 9b) represents the observed and modelled tidal ellipses from runbc for the M2 (M4) component. In most areas, the modelled tidal currents magnitudes are found to be weaker by up to 15 % with respect to the observations, but within the assumed error bounds. Also, the model seems to slightly underestimate semi-major axes almost everywhere, the RMS difference between observed and modelled tidal semi-major axes being 4.6 cm s −1 , again arguably within error bounds. The main differences between model and observations are found in the eastern (coastal) part of the domain, where islands channel and intensify tidal currents. It is, however, known that radar data quality is affected near islands and in the vicinity of the coast (Muller et al., 2009). In addition, Le Boyer et al. (2009) show that in that area the interpolated OMA outputs tend to significantly depart from the raw data. Therefore, more detailed or ad hoc studies will be required in order to better identify the source of the discrepancies whenever applications will require good estimates of tidal currents.

Tidal fronts
The occurrence of thermal fronts due to tidal mixing is a common feature on the European continental shelf. Thermal fronts appear in summer when tidal mixing prevents surface stratification induced by solar insolation from setting up. They separate tidally mixed waters from stratified waters. Following the study of Holt et al. (2008), we estimate the capability of the model in reproducing the location of tidal mixing fronts by comparisons with mean frontal positions as estimated from climatological ICES data (http://www.ices.dk/). We perform a similar analysis over 2008: data consist in all ICES temperature profiles available in June-August period, www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 and mean tidal fronts are deduced from the observed surface to bottom difference defined by T = 0.5 • C. The same approach is applied to model temperature profiles. Results are presented in Fig. 10. Over the shelf, SST fronts are due to the vertical turbulent mixing generated by tidal stresses at the seabed; the model shows reasonable agreement with observations in the main tidal front positions over the shelf, particularly for the Ushant front (northwest coast of Brittany), at the entrance of the English Channel, in the northern Irish Sea, off North Scotland and south of Doggar Bank (along ∼ 54 • N). MODIS SST data have also been used to perform qualitative comparisons of the SST fronts (Fig. 11). The modelled SST shows very good agreement with observations in the main tidal front positions. Over the shelf, the model presents larger stratified regions in the Irish Sea and colder SST in the Ushant front (off West Brittany) and at the entrance of the English Channel. Elsewhere SST fronts are relatively well represented both in temperature and in extension. On the Armorican shelf break, internal tides induce strong vertical mixing, which prevents the formation of seasonal stratification and results in a cold SST tongue along the shelf break from 11 • W to 5 • W. The cold intrusion is also reproduced in the model: temperature is underestimated by 0.5 • C and the southeastward extension is relatively well represented.
The relatively good accuracy of the model in reproducing these tidal mixing fronts is mainly related to the good representation of the bottom stress which affects the propagation of barotropic tidal waves,the k-epsilon two-equation model combined with the closure model of Canuto et al. (2001), the limitation of the dissipation under stable stratification and the tracers' advection (Holt et al., 2008). The model presents weaker stratification over Doggar Bank and to the west of England, and summer stratification is overestimated in the southern Irish Sea. Differences with observed mean frontal positions are also seen southwest of Denmark. Salinity stratification is important in this region (two major rivers including the Elbe river) and may impact frontal position (Holt et al., 2008).

Surges
For verification purposes, we diagnose the model's ability to reproduce surges by comparing residual elevations (i.e. detided sea levels) from the model with in situ measurements at time scales shorter than 10 days. At these time scales, the simulated coastal sea level is the signature of the static and non-static response to atmospheric pressure, of the response to wind forcing and of tidal interactions (see, for instance, Carrère and Lyard, 2003). The in situ sea level dataset consists of 115 tide gauges located along the European coasts (see red dots on Fig. 1). The harmonic analysis tool from Puertos del Estado is used to compute the tidal component from tide gauge elevations and to predict the tidal elevations; this software was developed to be used within the ENSURF Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ ensemble system (Pérez et al., 2012) and is based on Foreman (1977) harmonic analysis software. Comparisons between residual elevations from the model and the observations are shown in Fig. 12 for the Bay of Biscay and the English Channel. The elevations due to the static response of the ocean to atmospheric pressure (the so-called inverse barometer or IB effect) alone have been included in the comparison. The model shows very good agreement with data: 90 % (67 %), of RMS differences are smaller than 5 cm ( 3 cm), and 80 % of correlations are greater than 0.9 (results from detailed comparisons between observed and modeled values; results only synthesized on Fig. 12) and it performs better than the IB response (54 % with RMS difference smaller than 5 cm). It is most effective at high latitudes, which are more energetic areas, and where the ocean dynamical space and time scales are smaller, while at midlatitudes the model response to atmospheric forcing is closer to the IB approximation (not shown).

Surface tidal currents over the slope
Surface tidal currents from runbc are compared to observed currents at buoy stations along the northern Iberian coast for the M2 constituent (Fig. 9c). The model tends to overestimate the current amplitudes almost everywhere. The largest differences are found at Villano-Sisargas: the modelled ellipse azimuth is remarkably consistent with the observations, but the modelled semi-major axis is more than twice as large as the observed one. García-Lafuente et al. (2006) also found the M2 barotropic tidal currents were much greater than expected in this area; they attributed it to internal tides of considerable amplitudes. Pichon and Corréard (2006) showed that the northwest Spanish continental slope is indeed an area of internal tide generation. In the vicinity of Villano-Sisargas, differences of the M2 surface current amplitudes between the run with stratification (runbc) and the run without stratification (runbt) present high values due to internal tides (Fig. 9c). As internal tide generation is very sensitive to the bathymetry gradient and as the bottom topography and the coastline geometry are complex in the Cape Finisterre area, we expect errors in the slope position in the model bathymetry to explain the differences observed at this station. Without any better estimate of the local bathymetry, we cannot further test this www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 hypothesis. At Cabo Silleiro, tidal ellipse parameters are not correctly modelled. The model ellipse inclination is parallel to the isobaths, which is consistent with the study of Visser et al. (1994) in well-mixed conditions. Visser et al. (1994) also showed that in a region of freshwater influence, the stratification significantly influences the cross-shore component that can reach 40 % of along-shore component; this is consistent with the observed ellipse orientation. At Donostia and Matxitxako, modelled ellipses are also oriented along the slope, while observed ones are oriented across the slope. Indeed, the vertical profile of the observed M2 semi-major axis shows a clear stratification in the surface layers which is not represented in the model (Fig. 13). Below this layer (upper 30 m) the modelled profile is in good agreement with observations for the two current profiles. As discussed in Sect. 4.3.3, the surface layer's salinity suffers from an over-or underestimation in the model with respect to observations (Fig. 14). It should be stressed that model-data comparisons at single locations are very challenging in the sense that the modelled field can vary a lot from one grid point to the other, especially in the slope area. Indeed, a small local error in the bathymetry (depth and slope) can generate large displacements of the circulation patterns. Therefore, such comparisons are not necessarily representative of the overall model performance. Here, we presented them to illustrate some physical processes at work in order to explain the modeldata discrepancies.

Internal tides
We use 10 yr of TOPEX/Poseidon data combined with 7 yr of Jason-1 altimeter data spanning the period from September 1992 to January 2009 to obtain estimates of the internal tide signature in sea surface elevation. Altimetric data are not used for the validation of barotropic tides, as barotropic tides are compared to the FES2004 solution which assimilates altimetric data. The fraction of internal tides which is phase-locked with astronomical potential has sufficient coherence in both space and time for its surface signal to be detected in altimetric multi-year time series (Ray and Mitchum, 1997). We use an along-track harmonic analysis product provided by CTOH/LEGOS (F. Lyard, personal communication, 2010) to estimate the temporally coherent M2 component. We work with TOPEX/Poseidon/Jason1 track 137 that crosses the Armorican shelf slope at nearly a right angle ( Fig. 1) as we expect the M2 internal tides to propagate almost in the same direction (internal tides propagate perpendicular to the shelf break, Pairaud et al., 2010). To separate the barotropic tides from the time-coherent internal tide signal, we then spatially filter the M2 real and imaginary parts so that we remove the barotropic large-scale signal (Ray and Mitchum, 1997). According to Pairaud et al. (2010), the wavelengths of baroclinic mode 1 and mode 2 are respectively 141 km and 75 km in the abyssal plain. Note that these values have been obtained using data collected during the period September-October 1994 and may vary with the stratification. To take into account the angle between the altimeter Ocean Sci., 9, 745-771, 2013 www.ocean-sci.net/9/745/2013/ Fig. 12. Taylor diagram for comparisons between tide gauge data and the model residual elevations (blue) and the inverse barometer (red). Results are computed for periods below 10 days and for data located in the Bay of Biscay and in the English Channel. The radial coordinate indicates the standard deviation normalized by the data standard deviation; the concentric circles represent the RMS differences for time series normalized by the data standard deviation (the observations are indicated with the black point), and the angular coordinate represents the correlation with observations. track and the internal tide direction of propagation we use a spatial cut-off slightly larger than the wavelength. We also filter at shorter scales to reduce the altimeter noise. Finally, along-track tidal components are filtered between 140 km and 210 km for mode 1 and between 60 km and 110 km for mode 2. The "residual tides" are obtained by removing the filtered signal from the harmonic components. Figure 15 presents the along-track amplitude and phase as computed from the altimetric data and from the model for mode 1 and for mode 2. It also presents the along-track bathymetry. Indications attesting that the small oscillations are internal tides are as follows: (1) oscillations are generated in the Armorican shelf slope region (45.8 • N-46.0 • N); (2) south of the shelf break the phase increases constantly in the Bay of Biscay, which indicates that the wave propagates away from the shelf; and (3) for mode 1 and mode 2 we find wavelengths of approximately (λ 1 , λ 2 ) = (197 km, 76 km) in the altimetric data and the model; these values are slightly larger but in good agreement with theoretical values. The agreement between the altimetric and model estimates is remarkable: the misfits lie between 1 and 2 cm for the to- tal M2 amplitude and the model presents similar wavelengths than altimetric data. Interestingly, both the model and data indicate that the first mode does not have the dominant surface signature in the Bay of Biscay. Minor discrepancies in the baroclinic amplitudes and phases are, however, observed. For mode 1 (mode 2) there is a southward (northward) shift of the model sea surface signal. For the two modes, model amplitudes are slightly overestimated, especially near the generation site for mode 2. As previously mentioned, errors in the model bathymetry (here for the Armorican shelf slope representation) are expected to lead to uncertainties in the tidal representation, especially close to the generation sites of internal tides.

Local stratification off the Basque Country
The Donostia and Matxitxako deep-sea buoys deployed by AZTI-Tecnalía over the slope off the Basque Country provide hourly measurements of atmospheric parameters and ocean temperature and salinity from 10 to 200 m (Rubio et al., 2013;Fig. 1, Table 4). They represent very valuable data as they are the only available time series of temperature and salinity profiles over such a long period. The time evolutions of temperature and salinity as observed and modelled at the Matxitxako station are presented in Fig. 14 profiles show clear seasonal stratification. The model reproduces reasonably well the rapid heating and freshening of spring 2008. The cooling and mixed-layer deepening during winter are also correctly represented. However, the model underestimates the MLD equivalent in winter (The MLD equivalent is computed as the depth at which the density changes from the density at 10 m by 0.02 kg m −3 ; results not shown). In summer 2007 we also notice that the model salinity is too low in the upper layers. This is due to the initialization from PSY2V3, which is too fresh in this region. On the contrary, the model is slightly saltier in the surface layers during summer 2008. However, rapid freshening pulses during September 2008 are well represented in the model. Surface salinity uncertainties reflect mainly those on the precipitation fields; the latter are unknown but are expected to be larger in coastal areas such as northern Iberia because of orographic effects. Furthermore, Ferrer et al. (2009) have shown that fresh water discharge from several Spanish rivers can im-pact the offshore salinity conditions under specific oceanometeorological conditions. During the winter of 2007-2008, the ocean becomes warmer and saltier over the whole water column. This phenomenon is more marked in the model and affects the 350 m upper layer: it corresponds to a warm extension of the Iberian Poleward Current, also called the "Navidad event" (Pingree and Le Cann, 1992). This pattern is discussed later in Sect. 4.3.5.

Residual currents
Data used for comparisons consist of residual currents from the surface buoys of Puertos del Estado (Table 4).  Poleward Current along the northern Spanish coast. In summer, some current pulses are also observed in the data from mid-June to late August 2008, and are well reproduced in the model.

A focus on the circulation in the southern Bay of Biscay
We now focus on the ability of the model in reproducing a typical dynamic feature of the southern Bay of Biscay: the warm extension of the Iberian Poleward Current along the northern Spanish coast during winter. In this section our objective is not consistency assessment; we analyse several variables or quantities that are relevant for the Navidad occurrence (e.g. temperature and salinity anomalies in the upper layers, poleward along-slope currents) and show how the model can be used as a complement to available www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 observations. The analysis of the modelled SST and sea surface salinity (SSS) in the Bay of Biscay shows a warm and salty current flowing poleward along the western and the northern Iberian Peninsula and reaching the southwestern French coasts during winter of 2007-2008 (Fig. 17, SSS not shown; see also wintertime salinity anomalies on Fig. 14). The extension of the Iberian Poleward Current at or near surface along the northern Iberian coast is a recurrent feature of the autumn/winter circulation of the Bay of Biscay (Frouin et al., 1990;Garcia-Soto et al., 2002;Le Cann and Serpette; and has been referred to as the "Navidad event" (Pingree and le Cann, 1992). The extension of the IPC and the associated warming are, however, very variable from one year to the other; furthermore, it also undergoes some variability at timescales of a few days, as discussed below. Its main signature is a surface warming along the coast as it carries poleward warm water masses from the western Iberian region. Satellite observations clearly show the warm intrusion of the Iberian Poleward Current extending poleward as far as the Cantabrian shelf slope north of Spain on 20 January 2008 (Fig. 17). The intrusion is represented by a similar spatial extent in IBI with its leading edge reaching 3.5 • W at the same date. Modelled SST is slightly warmer than observed SST within the Iberian Poleward Current and in surrounding waters; moreover, the warm tongue appears to be narrower northwest of Spain. At the surface, the current thermal signature is about 1 • C warmer than adjacent waters in the data and in the model. It decreases poleward, varying from more than 14.5 • C off northwest Spain to about 13.5 • C at 3 • W. This corresponds to a drop of 1 • C over a ∼ 450 km distance, which is in good agreement with the poleward temperature gradient obtained by Frouin et al. (1990) using satellite SST data. The current narrows from 40 km off northwest Spain to less than 30 km at 3 • W. Vertical across-slope sections of thermohaline fields along the warm current (not shown) indicate the core of warm and salty water flows over the upper shelf slope, and extend from the surface to 350 m. The temperature field structure is illustrated with the map of 200 m temperature (Fig. 17).
Comparisons between modelled and observed along-slope residual currents (Fig. 16) indicate that the time variability of currents in the model is in reasonable agreement with data. They also suggest that the model tends to overestimate the IPC surface currents at the Puertos del Estado buoys (e.g. Cabo Silleiro and Estaca de Bares). At Bilbao (Fig. 16) the eastward pulses in February and March are not observed; this suggests that the IPC may be penetrating too far into the Bay of Biscay in the model. The analysis of other Navidad events (i.e. over other years) should be necessary to determine whether the model overestimation of the zonal extension of the IPC is a systematic bias in this configuration or if it occurred specifically for these events of February-March 2008.
The warm current transport has been estimated along the slope for sections perpendicular to the flow (see definition in the caption of Fig. 18). Considering typical temperature and salinity values along the current in the model, the transport has been computed for water masses warmer than 13.0 • C and saltier than 35.7 psu. Note that these values do not take into account the spatial gradient of temperature and salinity in the flow. The transport time evolution is shown in between 13 and 24 cm s −1 , which agree with typical values published (Frouin et al., 1990;Garcia-Soto et al., 2002). From Fig. 18a, we also notice the warm current transport is highly variable at periods of a few days and may vary from no transport to more than 2 Sv over a short period. Such variability at daily timescales has been described by Herbert et al. (2011) during the winter of 2004 from in situ observations and a numerical simulation. These authors suggest that the high-frequency variations are due to local wind forcing. When it is established, the current varies almost simultaneously at each section in response to forcing. The transport time series present a lag of 1 day between 8 • W and 6 • W (6 • W and 3.8 • W) and correlation of 0.75 (0.79) between the two transports. Driving mechanisms of the poleward current are often attributed to the large-scale meridional pressure gradient in the North Atlantic and to local wind stress (Frouin et al., 1990;Le Cann and Serpette;. Figure 18b represents the meridional wind stress off western Iberia and of the zonal wind stress off northern Iberia averaged over boxes defined following Le Cann and Serpette (2009). Positive values of wind stress indicate downwelling-favourable conditions which are likely to force the poleward flow. West of Iberia, wind stress is upwelling favourable until October 2007; then it is northward until the end of the year and is expected to provide a driving mechanism of the poleward current. When the current has reached Cape Finisterre and turns eastward, its poleward flow is maintained by eastward winds.

Conclusions
A high-resolution simulation covering the Iberia-Biscay-Ireland (IBI) region was carried out over the period of July 2007-February 2009 within the MyOcean project (Cailleau et al., 2012). In this paper, we had three main objectives: (1) to give an overview of the model configuration used for the simulation; (2) to give a broad-brush account of one particular aspect of this work, namely consistency verification; this type of validation is conducted upstream of the implementation of the system before it is used for production and routinely validated; it is meant to guide model development in identifying gross deficiencies in the modelling of several key physical processes; and (3) to show that such a regional modelling system has potential as a complement to patchy observations (an integrated approach) to give information on nonobserved physical quantities and to provide links between observations by identifying broader-scale patterns and processes. Objectives (2) and (3) were pursued in the year 2008.
We then examined domain-wide consistency verification results in terms of barotropic tides, transports, sea surface temperature and stratification. We also focused with more detail on two dynamical subregions: the Celtic shelves and the Bay of Biscay slope and deep regions; there, the modeldata consistency was checked, not exhaustively, but for a few important variables and processes such as tidal currents, tidal fronts, internal tides and residual elevation. Although a number of relatively minor points have been raised in the paper and some of these points are not in the IBI system specifications (e.g. details of the MW outflow), no major discrepancy with respect to expectations was identified. However, the study highlighted a few grey areas requiring more specific attention and testing, as well as priorities for improvement.
For instance, the analysis of the tidal elevations and current fields suggests that local improvements in the bathymetry are essential to better represent the tidal propagation, the bottom mixing and the generation of internal tides. Such improvements are also essential in regions of downscaling to finer-grid models, as well as to facilitate point comparisons with observations, a small error in the bathymetry (depth and slope) being able to generate large displacements of the circulation patterns.
Another important forcing in coastal areas is freshwater discharge, mainly from rivers. The analysis of the baroclinic www.ocean-sci.net/9/745/2013/ Ocean Sci., 9, 745-771, 2013 tidal surface signal points out at least two regions in the Bay of Biscay where the surface stratification in the model was probably poorly represented (however, the lack of available subsurface hydrological data will not permit direct comparisons). In these regions and according to past studies (e.g. Ferrer et al., 2009), runoffs from the Galician and the Basque rias have a significant impact on the local stratification but are not taken into account in our model configuration. Tests are necessary to quantify the impact of local changes of the surface stratification at different timescales on the baroclinic tidal currents. More generally, it becomes evident that a large effort for a more realistic representation of continental water discharge including the runoffs from small rivers is needed at least in these areas. Following some of our consistency checks, daily runoffs (instead of monthly climatological ones) are now prescribed in the operational version of this IBI configuration. Other points follow. Further increases in horizontal and vertical resolution would allow for some types of applications to be better served, such as biological connectivity and surface currents, the latter being probably strongly influenced by wave-current interactions (Law-Chune et al., 2013). Further tuning of the bottom friction and its spatial dependence on the seabed type may improve the solution locally for tides and surges, with impacts on SST distribution in frontal regions. Uncertainties in atmospheric forcings are likely to drive model errors on surface circulation and on tracer distributions via wind-driven turbulence. Errors are expected to occur for high frequencies and small-scale patterns, especially close to the coasts where the effects of smallscale orography are poorly represented. Stochastic modelling would provide an interesting framework for the assessment of the impacts of atmospheric uncertainties (e.g. Quattrocchi et al., 2013). We also suggest that R&D studies should be carried out to accompany the development of atmospheric downscaling for operational systems.
The model-data comparisons presented in this paper did not include an error analysis; for this, we need reliable estimates of uncertainties on the data, which are generally unavailable. Such estimates are required, for instance, in order to better identify the source of the discrepancies between the simulations and the gridded HF radar data in the vicinity of the Ushant front. The representativeness of local measurements (such as T /S profiles from the EN3 database or current time series from the Puertos del Estado buoy network) is a general problem, linked to the problem of upscaling local information when comparing to increasingly highresolution simulations. We also suggest to complement the datasets used for the consistency checks: for example, both Ferrybox T /S data and satellite altimetry data would be a nice complement, both datasets being characterized by repetitive sampling.
In the last section of this paper, we investigated the representation of the winter slope current along the northern Spanish coast (Navidad event). As very few data were available, this test is an illustration of the potential of our modelling system as a complement to observations within an integrated approach. We examined diagnostics on observed (e.g. surface currents) and non-observed (e.g. along-slope transport) quantities that could serve as a basis to define indices and that provide links between observations. Lastly, the relatively high (2-3 km) horizontal resolution used clearly pushes the model into the sub-mesoscalepermitting regime over much of the domain and allows for a significant part of the internal wave spectrum to be resolved. This seems to be the natural way forward to anticipate future high-resolution satellite observing systems such as SWOT (Surface Water and Ocean Topography).