Interactive comment on “ A coupled wave-3 D hydrodynamics model of the Taranto Sea ( Italy ) : a multiple-nesting approach ” by M . G . Gaeta

The present manuscript describes the implementation of embedded multiscale modelling systems in the Taranto Sea to refine predictions of hydrodynamics components (significant wave height, sea levels and currents) in the Mar Grande area. This research study provides, in particular, further insights about capabilities of a multiplenesting method to improve numerical predictions at the local scale of harbor basin which requires generally improved temporal and spatial resolutions. Of particular interest is the availability of numerous in-situ measurements in the area considered. Nevertheless, strong differences are exhibited between observations and results from


Introduction
The ability to accurately represent hydrodynamics processes in nearshore waters is essential for today's operational applications for coastal planning, decision support and assessment, since the majority of human activities are concentrated within the coastal zone.The estimate of sediment transport rates, the design of coastal defence structures and harbors, the forecast of pollutant concentrations, the regularization of ship routing and the knowledge of tide and wave energy potential at specific locations are all among the applications requiring a complete reproduction of the wave-and currentinduced processes.
High-resolution wave and hydrodynamics modeling offers an extensive range of capabilities regarding simulated conditions, works and practices, as well as providing a wide array of data regarding nearshore hydrodynamics.
Nowadays, predictive operational oceanography takes into account regional, subregional and shelf-coastal scales, based on coupled models of wave, current and baroclinic (sea temperature and salinity) dynamics, including interactions with the atmosphere by means of empirical bulk formulae (Guillou et al., 2013).However, since the above-mentioned oceanographic operational systems do not reach spatial resolutions lower than a few hundred meters, they are not able to take into account typical features of a coastal engineering scale: (i) nonlinear processes of wave propagation, (ii) interactions between offshore and coastal structures and (iii) boundary conditions induced by coastal processes, like river discharges (Saìnchez-Arcilla et al., 2014).The combined effect of waves and currents in shallow waters was furthermore the subject of numerous studies since Longuet-Higgins and Stewart (1964) introduced the term of radiation stress in hydrodynamics equations, a net momentum flux induced by spatial variation of the wave action spectra responsible for phenomena like wave setup, wave set-down and longshore currents.
Several authors have demonstrated that the coupling of wave and surge, and tide (Holthuijsen, 2007;Roland et al., 2009;Wolf, 2009) and ocean currents (Hersbach and Bidlot, 2008;Benetazzo et al., 2013;Brando et al., 2015), is a key element influencing the accuracy of nearshore dynamics predictions.The most relevant aspects of these interactions define the influence of water-depth variation in modulating wave energy during propagation and breaking, and of tidal currents in steepening incident waves.Hence, correct wave forecasting in coastal waters is only possible through the implementation of a multiple-nesting methodology, based on the development of (i) numerical models at high spatial and temporal resolutions, (ii) a downscaling technique and (iii) two-way coupling between hydrodynamics and spectral models.
Numerical models based on unstructured grids have the advantage of describing complicated bathymetry and ir-regular boundaries in shallow water areas more accurately and have recently been implemented to combine large-scale oceanic (regional) and small-scale coastal (local) dynamics (Lane et al., 2009;Liu and Xie, 2009;Xing et al., 2011;Guillou and Chapalain, 2012;Ferrarin et al., 2013), showing a significant improvement in predicting the wave field.The availability of field measurements is essential in the above context in order to calibrate and validate the numerical models, thus enhancing their predictive capabilities (Scroccaro et al., 2004;De Serio et al., 2007).
The objective of the present paper is the description of an operational strategy for the development of a multiscale modeling system, based on a multiple-nesting approach and open-source numerical models, and its implementation -at the highest resolution -of the dynamics of the Mar Grande at Taranto (southern Italy).Section 2 presents the methodology of the proposed multiple-nesting approach, describing the four levels of downscaling and the features of the adopted numerical models.The case study of the Taranto Sea is described in Sect.3, together with the available field measurements, used to calibrate and validate the adopted numerical models, thus confirming the validity of the proposed approach.In Sect.4.1, simulation results of 2-D wave-current interactions for the entire Gulf of Taranto are shown to provide the offshore boundary conditions in terms of wave forcings for the small-scale 3-D model of the Taranto Sea.The computed current pattern developing in Mar Grande and the evaluation of flux exchanges with the open sea are discussed in Sect.4.2, in comparison to the collected field data.

Multiple-nesting approach and model description
The proposed multiple-nesting approach is represented in Fig. 1 where the computational domains and bathymetries of the different models composing the numerical chain are shown.The methodology is based on four levels of downscaling as described in the following.
a.The large-scale system for the entire Mediterranean Basin (MFS, Mediterranean Forecasting System, e.g., Pinardi et al., 2003;Pinardi and Coppini, 2010), the first model of the chain (Fig. 1a).Madec et al., 1998) finite-difference code with a horizontal resolution of 1/16 • (7 km approximately) and 72 vertical σ −levels.The forecasting system performs a data assimilation system based on the 3DVAR scheme developed by Dobricic and Pinardi (2008).
b.The Southern Adriatic Northern Ionian coastal Forecasting System (SANIFS, in Fig. 1b), a coastal-ocean operational system providing short-term forecasts and built on the 3-D unstructured-grid finite-element hydrodynamic SHYFEM model (Umgiesser et al., 2004;Ferrarin et al., 2013) c.The coupled wave-2-D hydrodynamics model of the Gulf of Taranto (Fig. 1c) is based on the respective modules of the TELEMAC suite (Benoit et al., 1996).In its setup, waves are propagated from the Crotone buoy using the spectral module TOMAWAC (henceforth denoted as TOM; Benoit et al., 1996) and wave-induced currents are reproduced using the 2-D hydrodynamics module TELEMAC2D (henceforth denoted as TEL2D; Hervouet, 2007).The model is forced by SANIFS results (Fig. 1b) as initial/boundary conditions at both the offshore boundary and the free surface.
d.The coupled wave-3-D hydrodynamics model of the Taranto Sea (Fig. 1d) is also based on the respective modules of the TELEMAC suite (Benoit et al., 1996).For the spectral module, the imposed waves at the offshore boundary are extracted from the TOM results of the 2-D coupled model of the Gulf of Taranto (Fig. 1c); the hydrodynamics module TELEMAC3D (henceforth denoted as TEL3D, Hervouet, 2007) is driven by predictions extracted from SANIFS simulations (Fig. 1b) and COSMO-ME outputs (for 3-D offshore boundary and surface conditions, as well as the 3-D initial conditions).
Table 1 shows the numerical features of the proposed nesting approach: the simulated period, the mesh-size range and the spin-up time are reported for the different scale models.The spin-up time, i.e., the initial runtime needed in order for the models to reach relative stability regarding hydrodynamics, decreases with the scale of the implemented computational domain.
Table 1.Setup of the models used in the present multiplenesting approach: simulated period, mesh-size and spin-up time for MFS (Fig. 1a), SANIFS (Fig. 1b), TEL2D+TOM (Fig. 1c) and TEL3D+TOM (Fig. 1d) models.The development of the very high-resolution models (Fig. 1c and d) of the proposed downscaling chain is the main object of the present study, as the overall operational strategy.The TELEMAC suite is distributed under a General Public License (GPL) and is available at TELEMAC (2015).It is a finite-element-based solver for shallow water flows, wind wave propagation, ground water flows, tracer transport, sediment transport and morphodynamics.The parallel version of the last released suite was used on GALILEO, a CINECA machine characterized by Linux Infiniband Cluster architecture with 8256 nodes, operating under the EU Innovative Super Computing Research Allocation (ISCRA) program.

Models
In the proposed approach, the wave and 2-D/3-D hydrodynamics modules of the TELEMAC system are implemented in order to propagate offshore waves and currents and reproduce nearshore dynamics.
TOM is a third-generation spectral wave model, solving a simplified equation for the spectro-angular density of wave action (Booij et al., 1999) by means of a finite-elementtype method, in order to describe wave propagation towards coastal areas.
TEL2D simulates 2-D hydrodynamics by solving the depth-averaged flow equations proposed by de Saint-Venant (1871) in order to study free surface hydraulics and tracer transport for both transient and permanent conditions.The space is discretized by a series of Delauney triangular unstructured elements.Source terms in dynamic equations are modeled to represent the Coriolis force, bottom friction and wind action; active (temperature and salinity) and passive tracers can be also reproduced, as source or sink terms.
Coupled wave-2-D hydrodynamics modeling was discussed and validated in Samaras et al. (2016), based on the intercomparison of the TELEMAC suite and the well-known commercial software MIKE21 by the DHI © group.
TEL3D (Hervouet, 2007) solves the 3-D Navier-Stokes equations, with the option of the non-hydrostatic pressure hypothesis, and includes (i) the use of a finite element unstructured grid, which allows selective refinement of the mesh at key locations in the domain and boundary fitting method for vertical discretization; (ii) the transport-diffusion equations of intrinsic quantities (temperature, salinity, concentration), in order to reproduce 3-D hydrodynamics including the transport of active and passive tracers; (iii) a wide range of options for vertical turbulence modeling.The numerical solver used in the model is based on a fractional step technique in which the governing equations are split into fractional steps and treated using appropriate algorithms for the advection and diffusion of flow variables.The advection of velocities and water elevations is done with the semi-implicit Streamline-Upwind Petrov-Galerkin scheme, and the conjugate gradient method is used to solve the diffusive terms.
The modules TEL2D and TEL3D can be directly coupled (two-way coupling) on the same computational mesh as the spectral module TOM in order to reproduce the dynamics of wave-driven currents (while keeping in mind the extent of applicability of wave-averaged models; Bellotti et al., 2003); the gradients of the radiation stress induced by waves are computed using the theory of Longuet-Higgins and Steward (1964) as part of the hydrodynamics equations.The TELEMAC modules transfers the updated values of current velocities and water depths to TOM, while TOM solves the wave action density conservation equation, and returns the updated values of the wave driving forces acting on the current to the hydrodynamics modules (Hervouet, 2007).The effects of an enhanced bottom fiction coefficient due to wavecurrent interaction, as also experienced by Grant and Madsen (1979), are not included in the present release of the code.
In order to implement the proposed multiple-nesting approach, the authors properly modified each of the aforementioned modules in order for them to be able to read space-/time-varying conditions.Details of the forcings in the adopted models are reported in Table 2, describing the initial conditions (IC), offshore boundary conditions (OBC) and surface boundary conditions (SBC) of the modeled variables.
Blue Kenue ™ , a freeware tool developed by the National Research Council of Canada, is used to prepare the variabledensity triangular meshes of the study area and to visualize/process models' results.
The setup of TELEMAC applications consists of six sequential steps (to be applied for each module), namely mesh generation, adaptation of a steering file, extraction of initial/boundary conditions from larger scale models, scripting and compiling modules' subroutines in FORTRAN and postprocessing.
3 The case study

Description of the investigated area
The studied area (Fig. 1c), located in the Ionian Sea, is scientifically interesting due to its connection with the Adriatic Sea system, inducing the large-scale circulation structure known as the Taranto anti-cyclonic gyre, as well as other minor cyclonic structures near the shelf.In addition, the area is considered to be a vulnerable and sensitive area, affected by  chemical and biological pollutant industrial discharges and intense ship traffic (Di Leo et al., 2013).Accordingly, the accurate representation of its circulation regime is essential for the reconstruction of mixing and dispersion processes in the area as well (Mossa, 2006).
Figure 2 shows the study area of Taranto Sea, comprising Mar Grande, Mar Piccolo and the adjacent industrial area (see also Fig. 1d).The semi-enclosed basin of Mar Grande is delimited by Punta Rondinella, St.The bathymetry and shoreline used in the present work resulted from the digitization of nautical charts acquired from the Italian National Hydrographic Military Service.

Field data collection
Since the areas of Mar Grande and Mar Piccolo are of particular interest for their hydrodynamic regime and environmental fragility, they have been the set location of a series of short-and long-term monitoring activities over the years; the most recent ones are listed in the following, their characteristics also reported in Fig. 2. The instrument is able to provide, with a frequency of 30 min, information on (i) waves, in terms of significant height, peak and mean direction and peak period, (ii) currents and (iii) free surface elevations.The station also has instrumentation to record environmental variables of the atmosphere (i.e., wind speed and direction, air temperature and humidity, atmospheric pressure, net solar radiation) and sea (i.e., water temperature, pressure, conductibility, dissolved oxygen concentration).
-An RMN gauge (black square in Fig. 2) of the Italian Mareographic Network is located at 17  air temperature, atmospheric pressure and humidity and cloud cover at a 10 min interval.
In 2014, an oceanographic cruise campaign (Pinardi et al., 2016) named MREA2014 (Marine Rapid Environmental Assessment 2014) was performed in the area of the Gulf of Taranto during the period from 1 to 10 October, by using an Italian Navy Hydrographic Institute cruise ship to collect environmental data at three scales: large scale, shelf-coastal scale and coastal-harbor scale.At this latter one, which is of interest in this work, the following occurred.
-ADCP transects (green lines in Fig. 2) were carried out at five different sections (Punta Rondinella, western opening, Canale Navigabile and Porta Napoli and between the points of St. Paul Island and Capo S. Vito) to estimate water fluxes through the above-mentioned sections.
-Three groups of five drifters (purple stars in Fig. 2) were released in the period from 5 October 2014 at 10:40 to 6 October 2014 at 13:00, in order to estimate the surface currents (up to −2 m of water depth, Quattrocchi et al., 2016).
-31 CTD profiles (blue circles in Fig. 2) were collected on 5 October 2014 to measure water temperature and salinity profiles.
The drifter trajectories and ADCP mean velocities will be compared with the numerical results, and the reproduced patterns of surface velocities, including secondary circulation inside the basin, will be discussed in Sect. 4.  The Gulf of Taranto (Fig. 1c) was discretized using a variable density unstructured mesh in Blue Kenue ™ .Close to the shoreline, the reproduced mesh does not have a detailed resolution since the purpose of the 2-D model is mainly to propagate waves from offshore waters and to extract the wave field at the offshore boundary of the 3-D model; hence, harbor infrastructures, breakwaters and small islands are not represented in the area discretization.TEL2D and TOM are implemented on a domain extending the whole Gulf of Taranto from 16 • 28 to 18 • 19 E and 39 • 23 to 40 • 31 N.The unstructured mesh comprises 91 750 nodes and 179 000 finite elements with a size of 1 km offshore to 50 m close to the −10 m iso-depth.All the 2-D runs are carried out for the period 1-10 October 2014.
TOM is driven by the wave components extracted every 30 min from the offshore buoy located in Crotone (black circle in the right panel of Fig. 6).The wave spectra are discretized to 16 directions, with a bandwidth of 24 • .Nonlin-ear energy conservative processes are included in the model setup according to the authors' experience reported in Samaras et al. (2016).The model of Hasselmann et al. (1973) is used for bottom friction dissipation using a constant friction coefficient.The wind-generated spectral energy source is modeled using the WAM-cycle 4 model (Komen et al., 1984(Komen et al., , 1994;;Janssen, 1991) including the whitecapping dissipation term.The atmospheric forcing datasets from COSMO-ME (based on the COSMO model and providing outputs every 3 h; Baldauf et al., 2011, among others) are adopted in the TELEMAC models, properly extrapolated to reconstruct the wind field in the zone between the last wet mesh nodes of the high-resolution model and the actual coastline following Kara et al. (2007).An interpolation procedure is used to create the oceanic fields over the TELEMAC grid, as described by De Dominicis et al. ( 2013), and space-/time-varying wind conditions are imposed in the models.
A model calibration is performed, varying the values of the wind drag coefficient, C D , and the bottom friction included in the model in order to achieve a good agreement with the measurements in terms of depth-average velocities, free surface elevation and wave height.The final setup of TEL2D includes a wind drag coefficient, C D , equal to 3.5 × 10 −3 in order to compensate the underestimation of the actual wind by the COSMO-ME fields.The bottom friction is computed according to Chezy's law, with a Gauckler-Strickler coefficient equal to 40 m 1/3 s −1 ; a constant Coriolis coefficient f is set equal to 1.101 × 10 −4 N m −1 .The Smagorinsky model (Smagorinsky, 1963;Hervouet, 2007) is chosen to represent the horizontal turbulent dissipation, with a velocity diffusivity coefficient equal to 0.1 m 2 s −1 .The final values for these coefficients were also used in the TOM module.Wave-driven radiation stresses computed in TOM are fed into the hydrodynamics module TEL2D, directly coupled in order to reproduce wave-driven currents.
Initial conditions of water elevations and depth-averaged currents are extracted from the SANIFS model and interpolated over the 2-D TELEMAC grid, imposing the no-slip conditions on the coastline.TEL2D is driven by the freesurface elevations and the depth-averaged currents (extracted hourly by SANIFS) and by wind forcings varying in time and space.A spin-up time of 2 days (1-2 October 2014) was necessary for the coupled model to reach stable conditions for currents.

Coupled wave-3-D hydrodynamics model of the Taranto Sea
The 3-D model of the Taranto Sea is implemented on a domain extending from 17 • 7 to 17 • 12 E and from 40 • 27 to 40 • 29 N (Fig. 3, top panel).The study area includes Mar Grande (maximum depth is −50 m), with all the islands and the coastal defence structures properly represented in the constructed mesh, including the harbors of the city of Taranto (Fig. 3, bottom panels) and the western open sea, in front of  the industrial area of the city, and Mar Piccolo, with the two small basins.
The horizontal unstructured mesh is constructed adopting the following mesh-size rules: an edge growing ratio of 1.2 : 100 m at the offshore boundary, where each node corresponds to a node of the SANIFS mesh, 50 m inside the liquid domain and 10 m at the coastline, harbor structures and islands.3-D space is represented as a series of 2-D triangular unstructured meshes with 55 399 nodes and 105 751 elements, at defined levels between the bottom and the free surface.The implemented vertical discretization follows a 17 z level scheme, with the setup adopted in SANIFS, where the thickness of each layer is equal to 2 m.
In TEL3D, three fractional steps are implemented to solve the governing Navier-Stokes equations; in the present study, horizontal turbulence is modeled using the subgrid Smagorinsky scheme (Smagorinsky, 1963;Hervouet, 2007); at the vertical scale, the mixing-length turbulence closure is applied following the classical formulation by Prandtl (1925).The values of bottom friction, wind drag and Coriolis coefficients used in TEL3D are equal to the ones estimated for the 2-D simulations of the Gulf of Taranto.Water temperature and salinity are modeled as active tracers (barocline circulation), following the mass-balance equation; their values at each time step and node are used to compute the water density variation according to the state equation (Hervouet, 2007).Waves resulting from the coupled TEL2D+TOM run of the entire Gulf of Taranto (see the following Sect.5.1) are imposed on the offshore boundary of the 3-D model.The subroutine, defining the imposed spectral energy conditions at the boundaries, was properly modified in order to reach each node of the offshore boundary and at each time step of the run the values of significant wave height, peak period and mean direction, extracted and linearly interpolated in time by the 2-D simulation results.
In analogy to the 2-D model, wind data from COSMO-ME are also imposed at the surface to include forcings that vary in time and space.
The spectral module TOM is coupled with the 3-D hydrodynamics model TEL3D, properly modified in order to account for time-/space-varying conditions (Table 2).Initial conditions of water elevations, velocities and tracers are extracted from SANIFS results on 3 October 2014 and are extrapolated and interpolated at each node of the domain and at each horizontal layer, adopting the procedure described in De Dominicis et al. (2013) to impose the data on the TELEMAC grid.A spin-up time of 2 days (3-4 October 2014) was necessary for the coupled model to reach stable conditions for both currents and tracers.
The three most significant waste-water discharges in the study area (De Pascalis et al., 2016) were also introduced/prescribed as boundary conditions in the model setup: the ILVA water pump withdrawal with yearly averaged discharge of about 27 m 3 s −1 from Mar Piccolo, and the two ILVA wastewater releases, northern Punta Rondinella and outside Mar Grande, with a total yearly averaged discharge of 26 m 3 s −1 .The contribution of submarine freshwater sources inside Mar Piccolo, which is not studied here, was neglected in the present study.Numerical results from the 2-D hydrodynamics module TEL2D coupled to the spectral module TOM are obtained in order to force the 3-D model of the Taranto Sea with the waves propagated from Crotone to the shelf-coastal area.In order to validate the results, the time series of free surface elevation and significant wave height are compared to the data collected at the Mar Grande station.The simulated period was characterized by small waves in Mar Grande until the day 5 October 2014 when wind increased and changed in direction.
Figure 4 shows a comparison at the Mar Grande and RMN stations of the wind velocity magnitude (top panels) and direction (bottom panels) between the measurements (red squares) and COSMO-ME values (blue lines) interpolated over the TELEMAC grid.The agreement is quite good, with the exception of the maximum event corresponding to the higher wave.This is attributed to the fact that a large-scale atmospheric model such as COSMO-ME does not properly reproduce local atmospheric events, mainly responsible for fast perturbation, such as the one occurring during the night between 6 and 7 October.
TEL2D managed to capture the tidal constituents well (Fig. 5), both in amplitude and frequency in comparison to the observations.Indeed, the harmonic analysis of the wa- ter elevation on an hourly basis was performed using the TAPPY tidal analysis package (Cera, 2011), and the results for the most important semidiurnal and diurnal constituents (M2, K1) were compared to observations.This analysis reports errors of 6.9 and 5.8 % for amplitude and phase of the K1 component, and errors of 3.7 and 2.2 % for amplitude and phase of the M2 component, showing very accurate results in comparison with other reference studies (Ferrarin et al., 2013;Federico et al., 2016).
The effects of tide-induced water elevations and currents on wave heights are analyzed in the area of the Gulf of Taranto, where waves measured at the Crotone buoy are propagated from the offshore boundary (Fig. 6).
The coupled TEL2D+TOM model predicts time series of wave height in comparison with measurements at the Mar Grande station satisfactorily, while the wave climate is roughly underestimated when using the standalone TOM run, as described in Fig. 7.The 2-D model results are in accordance with the field evidence reported in De Serio and Mossa (2016), who showed how tidal flow induces an increase in the wave height.Numerical results tend to slightly overestimate the wave height after the peak (where H s equal to 0.4 m) characterizing the storm between 6 and 7 October 2014; wind directions are shown to rapidly change during this short event, which the imposed atmospheric forcings from COSMO-ME does not reproduce.
Figure 8 (left panel) shows the color map of the significant wave height, H s , where x axis reports the time, y axis the relative distance on the offshore boundary.This output from the 2-D model is used as boundary conditions for the spectral module in the 3-D computation, as described in the following subsection.The prescribed signal variation on the boundary outline is extracted for the exemplary instance of 6 October  2014 at 14:00 (noted as a dotted line in the left panel) and is shown in the right panel.

Coupled wave-3-D hydrodynamics model of the Taranto Sea
The capability of the proposed multiple-nesting approach and of the implemented models is investigated on the reproduced dynamics of Mar Grande, focusing, in particular, on current trajectories and intensity.Figure 9 shows the color maps of the surface velocity field (top) and of the significant wave height (bottom) in the area of Mar Grande.In particular, the values represent the average occurring in the interval from 10:00 to 12:00 on 5 October 2014; the first two layers' depth averages (up to −2 m of water depth) of the 3-D velocity results were obtained in order to compare the numerical outputs with the average observed trajectories of each drifter group (black arrows in the top panel).The current pattern in the area where the drifters were released is adequately captured by the coupled model: the average objects' propagation is estimated to be directed towards around 270 • N, except for the drifter group D2, moving northwards, while the computed surface velocity is slightly underestimated in the area close to the release point of the drifter group D1. Figure 10 shows the depth-averaged velocities along four opening sections (green lines in Fig. 2; the ADCP data collected along the southern opening were neglected due to noisy signals for the severe climate conditions during the measurements) where the ADCP collected current profiles in the measurement interval from 07:20 to 09:30 on 5 October 2014.An excellent agreement between observations and the numerical results is observed for the mean flux between Mar Grande and Mar Piccolo (along Porta Napoli and Canale Navigabile sections) and the western open sea (along Punta Rondinella and western opening sections), with a mean velocity root mean square error value lower than 7 cm s −1 .The effect induced by the counterclockwise circulation at the surface in Mar Grande is also in accordance with the comparison at the ADCP transects as reported in Table 3: a negative mean flux exiting the Mar Grande basin is observed both in measurements and computations; discharges along the sections are also well reproduced by the model, along with the mean flux velocities and directions.The presented results also agree, as well as with the available measurements, with the findings by De Pascalis et al. (2016), that described the 2013-averaged fields of Mar Grande, while SANIFS outputs in Federico et al. (2016) reported opposite circulation and fluxes.
The time series of currents and the vertical profiles at the point of the Mar Grande station are presented in Figs.11 and 12, respectively: in both panels, the numerical results from the TEL3D+TOM run are compared with the observations and the SANIFS outputs.The overall agreement of computations is acceptable, with a mean error in the depth-averaged velocity magnitude equal to 24 and 33 % for TEL3D+TOM and SANIFS runs, respectively.Vertical profiles of velocity at the Mar Grande station are also well captured by the model: the daily-averaged values in the period from 5 to 7 October 2014 show a good agreement with observations, with decreased velocity BIAS values for TEL3D+TOM results, especially after 6 October.
Finally, the combined effects of (i) the coupling of hydrodynamics and wave dynamics; (ii) the inclusion of Mar Piccolo dynamics in the modeling system, (iii) the introduction of the main wastewater discharges; and (iv) the increase in spatial and temporal resolution have led to an accurate representation of the circulation regime inside Mar Grande by the proposed modeling approach, including the dynamics influencing the Mar Piccolo area.

Conclusions
The present work proposes a multiple-nesting methodology for the assessment of nearshore wave dynamics and hydrodynamics, based on a downscaling approach -from an oceanographic to a coastal engineering scale -and the implementation of open-source numerical models, customized for specific user purposes and accounting for wave-hydrodynamics processes.The approach, deemed to be useful for environmental planning and decision support systems for maritime hazards, was applied to scale down from the Mediterranean Forecasting System, MFS (regional system), to the Southern Adriatic Northern Ionian coastal Forecasting System, SAN-IFS (mesoscale shelf-coastal), to very high-resolution applications of the models comprised in the TELEMAC suite.
The aforementioned methodology was applied to the case study of the Taranto Sea, southern Italy, a vulnerable and sensitive area including Mar Grande and Mar Piccolo and the adjacent industrial area.Field measurements, available from a recent oceanographic campaign (MREA campaign, October 2014) and monitoring stations in Mar Grande, were used to calibrate and validate the TELEMAC models.Nesting regional and local systems based on the coupling of wave propagation and circulation modules have been used to propagate waves from the offshore boundary of the Gulf of Taranto to nearshore waters.2-D simulations of wavecurrent interactions were implemented through the two-way direct coupling of the spectral module TOMAWAC and hydrodynamics module TELEMAC2-D.Models were forced at the offshore boundary by (i) wave spectra from the Cro-tone buoy, (ii) free surface elevation and currents extracted from SANIFS and (iii) at the surface by wind data from the COSMO-ME model.The computed time series of sea levels and wave heights show a good agreement with measurements carried out in Mar Grande, revealing that the inclusion of tide-induced time-/space-varying water elevations and velocities in wave modeling is essential for an accurate representation of the wave field during storms.
A coupled wave-3-D hydrodynamics model of the Taranto Sea was applied to investigate the dynamics of the area of Mar Grande, focusing, in particular, on current trajectories and intensity.The consistent trends of the current field and circulation are well reproduced with the development of a small cyclonic pattern at the surface, in accordance with the observed trajectories of the released drifters.The fluxes of water masses between the open sea and Mar Grande and Mar Piccolo are well reproduced by the model too, showing a good agreement with the ADCP data.The intercomparison with SANIFS results shows an overall satisfactory agreement with the observations in terms of current profiles at Mar Grande station, and is deemed to provide useful insights on implemented model features.
The results highlight the capabilities of the proposed multiple-nesting approach to comply with the observed sea conditions at a coastal-harbor scale, considered to be essential for a wide range of coastal engineering and management purposes.Overall, the methodology is deemed to be of general interest for ocean and coastal modelers involved in (i) the development of procedures for nesting high-resolution models to general circulation models of regional/subregional scale (as well as their validation based on field data); (ii) the adaptation of coastal models to platforms of operational oceanography; (iii) the implementation of the above in coastal planning and design (coastal zones and harbors).

Data availability
The model and observational data used to carry out this research are free and available on request by writing to the authors.

Figure 1 .
Figure 1.Multiple-nesting approach from the Mediterranean Sea scale with (a) MFS and (b) SANIFS to a coastal-engineering-scale model performed with the TELEMAC suite (c) in 2-D for the Gulf of Taranto and (d) in 3-D for the Taranto Sea: maps show the adopted bathymetry.

Figure 2 .
Figure 2. Map of the Mar Grande study area and points/locations of available measurement datasets from the MREA14 campaign (CTD, ADCP and drifters), the Italian Mareographic Network (RMN) and the Mar Grande station (MEDA).
Peter and St. Paul Islands and Capo S. Vito and is directly connected to the Gulf of Taranto through one large opening at its south boundary and two smaller ones at the western boundary.Mar Grande hosts the commercial and military harbors of the city of Taranto, as well as a number of mussel farms close to Capo S. Vito.Mar Piccolo is divided into two basins, First Seno and Second Seno.It is characterized by the presence of a number of submarine fresh water springs, with a combined discharge of around 5 m 3 s −1 (De Pascalis et al., 2016), and by numerous mussel farms.It is connected to Mar Grande through two navigable openings, named Porta Napoli and Canale Navigabile.A submarine freshwater source with an average discharge of 17 m 3 s −1 is also present at the mouth of Porta Napoli (De Serio and Mossa, 2016).

-
The monitoring station in Mar Grande -MEDA -(De Serio and Mossa, 2016) is located at 17 • 12.9 E and 40 • 27.6 N, in the middle of the Mar Grande basin (red square in Fig. 2) and was installed in February 2014 under the Italian projects PON no.H51D11000050007 and RITMARE.The station consists of one upwardlooking acoustic Doppler current profiler (ADCP) placed on the sea bottom, at a depth of −22.65 m.

Figure 3 .
Figure 3. 3-D model of the Taranto Sea: the overall grid (top panel) and two detailed zones (bottom panels) close to the commercial harbor (a) and the military harbor (b) of Taranto; background images are from Google Earth, 2016, privately processed.

Figure 4 .
Figure 4. Wind characteristics at Mar Grande (left panels) and RMN (right panels) stations: velocity magnitude (top panels) and direction (bottom panels) from the interpolated COSMO-ME values on the TELEMAC grid (blue line) and from observations (red line and circle markers).

Figure 5 .
Figure 5.Time series of the tidal water elevation at the Mar Grande station: numerical results from the coupled TEL2D+TOM run (blue line) and observations (red line and circle markers).

Figure 6 .
Figure 6.Numerical results from the TEL2D+TOM run for the Gulf of Taranto: significant wave height from Crotone buoy (black circle in the right panel) imposed on the offshore boundary (left panel), and snapshot (right panel) of wave propagation on 6 October 2014 at 14:00 (color map of the significant height, H s , and arrows of the mean direction).

Figure 7 .
Figure 7. Significant wave height, H s , at the Mar Grande station in the period 5-7 October 2014: numerical results from the standalone TOM run (black dotted line), the coupled TEL2-D+TOM run (blue line) and observations (red line and circle markers).

Figure 8 .
Figure 8. Color map of the significant wave height, H s , imposed on the offshore boundary of the 3-D model in time (left panel) and prescribed signal variation on the boundary outline extracted for an exemplary instance (right panel): instance of 6 October 2014 14:00 noted as a dotted line in the left panel.

Figure 9 .
Figure 9. Numerical results from coupled TEL3D+TOM in the area of Mar Grande on 5 October 2014 (average over 10:00-12:00).(a) Surface velocity: color map for the velocity magnitude, white arrows for current directions and black arrows for drifter routes (average).(b) Wave distribution: color map for the significant wave height H s and black arrows for the mean direction.

Figure 10 .
Figure 10.Depth-averaged velocity along four opening sections in Mar Grande: TEL3D+TOM results (blue line) and ADCP measurements (red line and circle markers).

Figure 11 .
Figure 11.Depth-averaged velocity at the Mar Grande station: numerical results from the coupled TEL3D+TOM run (red line), SANIFS run (black line) and observations (red line and circle markers).

Figure 12 .
Figure 12.Daily-averaged vertical profile of current at the Mar Grande station and BIAS.Left panels: numerical results from the coupled TEL3D+TOM (blue line), SANIFS (black line) and observations (red line and circle markers); right panels: continuous lines show the velocity BIAS for the two models over the water depth.

Table 2 .
Forcings for the different variables modeled in the proposed multiple-nesting approach.
a Temperature and salinity in TEL3D.b Significant wave height, peak period and mean direction defining a JONSWAP spectrum in the module TOM.c Atmospheric pressure, air humidity and temperature, cloud cover and wind in TEL3D.

Table 3 .
Water mass exiting Mar Grande on 5 October 2014: comparison between numerical results (Num) from TEL3D+TOM and observations (Obs) from the four ADCP transects as reported in Fig.3.