Influence of Wave State and Sea Spray on the Roughness Length : Feedback on Medicanes

Occasionally, storms that share many features with tropical cyclones, including the presence of a quasi-circular “eye” a warm core and strong winds, are observed in the Mediterranean. Generally, they are known as Medicanes, or tropical-like cyclones (TLC). Due to the intense wind forcings and the consequent development of high wind waves, a large number of sea spray droplets—both from bubble bursting and spume tearing processes—are likely to be produced at the sea surface. In order to take into account this process, we implemented an additional Sea Spray Source Function (SSSF) in WRF-Chem, model version 3.6.1, using the GOCART (Goddard Chemistry Aerosol Radiation and Transport) aerosol sectional module. Traditionally, air-sea momentum fluxes are computed through the classical Charnock relation that does not consider the wave-state and sea spray effects on the sea surface roughness explicitly. In order to take into account these forcing, we implemented a more recent parameterization of the sea surface aerodynamic roughness within the WRF surface layer model, which may be applicable to both moderate and high wind conditions. The implemented SSSF and sea surface roughness parameterization have been tested using an operative model sequence based on COAWST (Coupled Ocean Atmosphere Wave Sediment Transport) and WRF-Chem. The third-generation wave model SWAN (Simulating Waves Nearshore), two-way coupled with the WRF atmospheric model in the COAWST framework, provided wave field parameters. Numerical simulations have been integrated with the WRF-Chem chemistry package, with the aim of calculating the sea spray generated by the waves and to include its effect in the Charnock roughness parametrization together with the sea state effect. A single case study is performed, considering the Medicane that affected south-eastern Italy on 26 September 2006. Since this Medicane is one of the most deeply analysed in literature, its investigation can easily shed some light on the feedbacks between sea spray and drag coefficients.


Introduction
In the last decades, satellite images have revealed the presence of intense vortices in the Mediterranean that show some visual, dynamic and thermodynamic similarities with tropical cyclones, Atmosphere 2018, 9, 301 2 of 17 like the presence of a central "eye" strong rotational winds, a warm core, weak vertical wind shear, a cold upper level trough in the early stages.These cyclones are generally known as Medicanes (acronym for Mediterranean hurricanes) [1] or Mediterranean tropical-like cyclones (TLC).
Different studies have analysed the properties of this category of cyclones [2][3][4].Baroclinic instability is important in the early phase of their lifetime, when they generally show extra-tropical features [1], while barotropic instability may play a relevant role in their mature stage [2].Their development depends critically on the intensity of sea surface fluxes and latent heat release associated to convection, like for tropical cyclones [5,6].In the low levels, the presence of a warmer sea surface temperature makes a larger amount of energy available for their development: however, this effect is case dependent [7,8].At upper levels, the presence of a potential vorticity anomaly appears to influence their evolution, in particular during the first phase [9][10][11].The conditions favourable for their development have been analysed from a climatological perspective by Cavicchia et al. [12] and Nastos et al. [13].The evolution of their characteristics, and in particular the occurrence of a tropical transition, can be analysed either using the Hart diagram [7,14,15] or by comparing the advective and diabatic terms in the surface pressure tendency equation [16].
Under strong wind conditions, like those occurring in Medicanes, the extreme wind forcing at the sea surface and the resulting development of high wind waves, with relative short periods and lengths, may produce a large quantity of sea spray droplets.In these conditions, tearing of the wave crests by the wind is an additional sea spray production mechanism other than the most common bubble bursting.The latter one occurs when bubbles, previously entrained in the water column by breaking waves, rise to the surface and burst, creating sea spray droplets with radius at formation (r 0 ) from submicron dimension to a maximum of less than 100 µm [17].The sea spray droplets generated via mechanical tearing are called "spume droplets" and their radius at formation is typically assumed in the range (20-500) µm with a peak centred at r 0 = 100 µm [17].A wind speed higher than approximately 7-11 m s −1 is typically necessary to generate spume droplets [18].
The effect of these spume droplets may be remarkable under several aspects.They change the basic state of the air-sea interface in the marine surface boundary layer (MSBL), affecting different physical processes that occur at the interface.Considering the combined action of atmospheric turbulence and gravitational settling, the largest sea spray droplets can be suspended in air only for a short period of time (few seconds), while the smallest droplets can stay longer (hours).Nevertheless, it appears that the largest spray droplets contribute mostly to the air-sea exchange of momentum and heat [19,20].Powell et al. [21] speculated that the increased foam coverage resulting from intensively breaking waves could progressively form a 'slip' surface at the air-sea interface that leads to the observed reduction of the sea drag coefficient at high wind speeds.This reduction has been found also in theoretical and numerical studies [22][23][24].
A common practice in numerical weather prediction models (NWP) is to calculate the air-sea momentum flux through the classical Charnock relation [25], which does not explicitly consider the wave-state and sea spray effects on the sea surface roughness.In order to take into account both effects, which may affect the roughness significantly, Liu et al. [26,27] and Wan et al. [28] used a parameterization of sea surface aerodynamic roughness length (z 0 )-based upon the work by Volkov [29] and Makin [24]-which may be applicable to both moderate and high wind conditions.
In the present paper, a single case study is considered in order to study the effect of both sea spray and wave-state on the sea surface roughness for the simulation of a tropical-like cyclone in the Mediterranean Sea.These two forcings cannot be neglected, taking into account the high wind speed reached during the event; in particular, we propose the inclusion of these two terms in the surface layer parameterization scheme of the WRF model, thus providing a more realistic simulation of the air-sea interaction processes.An operative model sequence based on the Coupled Ocean Atmosphere Wave Sediment Transport (COAWST) system [30] and the WRF-Chem model version 3.6.1 [31], has been implemented.First, the wave state conditions are provided to WRF via a wave model, SWAN [32], two-way coupled with WRF as part of COAWST, which includes the two-way interaction of air, ocean and wind waves [33].In a second stage, the effect of sea spray is also calculated, including a newly implemented Sea Spray Source Function (SSSF) in WRF-Chem [28].WRF-Chem is one-way coupled with COAWST: first, COAWST transfers the wave-state data to WRF-Chem; then, the GOCART [34] aerosol sectional module in WRF-Chem explicitly considers the emission and transport of sea spray aerosols via the implemented SSSF.Finally, the wave and sea spray data update the roughness length in the surface layer parameterization used in WRF-Chem.
The rest of the paper is organised as follows.In Section 2, "Materials and Methods," in addition to the model setup, we provide a short discussion of the synoptic conditions during the case study, together with the description of the mechanisms of sea spray generation under hurricane conditions and of the sea surface roughness parameterizations.Section 3 presents the results of the SWAN/WRF-Chem outputs, also considering the evolution of the cyclone along the track.Concluding remarks are in Section 4.

The Event: Synoptic Description
The tropical-like cyclone crossing south-eastern Italy on 26 September 2006 has been one of the most investigated Medicanes in the literature [6,35].Thus, an intensive study of the effects of waves and sea spray on the drag coefficients is possible, comparing the results of our integrated numerical approach with the deeply analysed features of this cyclone.
The track of the cyclone is shown in Moscatello et al. [6].First, the cyclone developed in the lee of the Atlas Mountains on 24 September 2006; next, the cyclone moved across the Strait of Sicily, where it intensified mainly due to strong latent heat fluxes from the warm sea surface on 25 September.Sensitivity experiments showed the importance of a warm SST for the intensification of this cyclone [36].After approaching the south-eastern tip of Sicily in the evening of 25 September, the cyclone started to follow an arc-like track, steered by the upper level flow, around the synoptic low, which remains nearly stationary over the Tyrrhenian coasts of central Italy (500 hPa, Figure 1a).
The mean sea-level pressure minimum strongly deepened in the early hours of 26 September, when the vortex started to cross the Ionian Sea near the east coast of Calabria.In this phase, intense convection developed, with high rainfall amounts recorded along the Ionian coast of Calabria: sensitivity experiments have demonstrated the key role of latent heat release due to intense convection in this phase [6].The crossing of the cyclone from the right to the left exit of an upper level jet stream (300 hPa, Figure 1b) at about 0600 UTC contributed to the intensification of convection [37].Thus, the interaction with the jet stream reinvigorated the cyclone, showing some similarities to the case study analysed by Reale and Atlas [2].
At about 0900 UTC, the cyclone reached the Ionian coast of the Apulia region, where the synoptic station at Galatina airport (G label in Figure 2) measured a wind gust of 78 knots [6] and a deepening of the pressure low by about 13 hPa in 90 min (not shown); a minimum sea-level-pressure of 986 hPa was recorded at Nardò station (belonging to a regional meso-network), near the Ionian coast of Apulia region (N label in Figure 2).Along the track in the Adriatic Sea, the cyclone showed the symmetric, deep warm core structure typical of a tropical cyclone: in this stage, both the surface fluxes and latent heat release appear important to self-sustain the vortex.Shallow convection was observed in this phase, which is a behaviour commonly reported in the mature stage of Medicanes [7,38].Finally, the cyclone made landfall in the northern part of Apulia region (Gargano promontory) at around 1700 UTC, when it was still intense (988 hPa).
The large-scale conditions, emerging from NCEP/NCAR reanalysis, show an environment favourable to the development of convection in the Ionian Sea: the high columnar precipitable water (Figure 1c) provide high humidity content in an extensive region; the steep lapse rate, due to warm temperatures in the low levels associated to a thermal ridge in the Ionian Sea and cold temperature in the upper levels, makes conditions favourable to potential instability, with Lifted Index locally less than −8 K (Figure 1d).Finally, the intrusion of Saharan dust in the upper level associated with a humidity content minimum [39] appears as a common feature for all these cyclones, whose role needs anyway still to be disentangled [11].The predictability of this cyclone was also analysed in a series of papers [37,40,41], which showed that the large-scale forcing (initial and boundary conditions), the numerical model used for the simulations and the choice of the microphysical scheme have all a significant impact on the cyclone track and its intensity.The use of a coupled atmosphere-wave model [42] should improve the simulation skill, through the use of a more refined sea surface temperature field, which also evolves realistically with time and through a more accurate representation of heat and moisture fluxes, consistent among the different modules of the numerical system.

Roughness Length over the Sea Surface
The momentum transfer between the atmosphere and the sea surface is generally described in terms of wind speed variation with height, a drag coefficient increasing with sea surface roughness (z 0 ) and the surface friction velocity (u * ).
The classical Charnock relation [25], which is valid for wind speed less than about 30 m s −1 , involves a quadratic relationship between the roughness length and surface friction velocity: where g is the gravitational acceleration and α is the Charnock parameter.The default option for the WRF model considers that the Charnok parameter is varied linearly from 0.011 to 0.018, between 10-m wind speed in the range [10,18] m s −1 and wave state effects are not explicitly taken into account.
To include both wave state and sea spray effects on the momentum roughness length, we have considered the parameterization of Liu et al. [26,27]: This is a combination of the resistance law of Makin [24], which reflects the impact of spray droplets on the airflow dynamics and the Volkov model [29], which depends on the wave age and is based on the observations of the Scientific Committee on Oceanic Research (SCOR) workgroup 101 [29].In Equation (2), the parameter β = c p /u * represents the wave age, which is calculated considering the ratio between the phase speed at the peak of the wave spectrum (c p ) and the surface friction velocity.The dependence on the sea spray is provided by the parameter ω = min(1, α cr /ku * ), which represents the correction on the logarithmic wind speed profile in the surface layer due to the presence of sea spray droplets (k ≈ 0.4 von Karman's constant).The critical value of the terminal fall velocity (α cr ) is set equal to 0.63 m s −1 according to Makin [24].Equation (2) may be utilized for both low-to-moderate and extreme wind conditions.
For both formulations described above, Equations ( 1) and (2), we have calculated the neutral drag coefficient as . The continuous black line in Figure 3 represents the drag as calculated by Charnock relation in Equation (1) with the Charnok parameter set to 0.018; the well-known unrealistic linear increase is apparent at the highest wind speed [43].The red, black dotted and blue lines in Figure 3 have been calculated with Equation (2), just considering three values of the wave age: β = 20, β = 27 and β = 50.Roughly speaking, small wave age values correspond to young sea states, while increased values signify fully developed sea states; if β = 27 Equation (2) reduces to the Charnock Equation (1) with a typical value of α = 0.018.In Figure 3, it can be seen that: (i) at wind values greater than the hurricane force, which correspond to friction velocity values greater than u * cr = 0.63/k = 1.54 m s −1 , the drag coefficient is levelled off and starts to decrease with a further increase in wind speed; (ii) in the considered wave age range, the mechanical drag at the sea surface decreases as the wave age increases.A possible explanation for the well-known decrease of the drag coefficient with the wind speed increase above the hurricane force (a critical threshold of about 33 m s −1 is typically assumed) is the development of a sea foam layer at the air-sea interface [21].In general, in the literature, the dependence of the drag coefficient on wave age is uncertain, because the slope of the fitted line can be influenced considerably by the narrow range of wave age [44].In addition, the presence of swell can generate complex interactions [23].Nevertheless, we believe that the behaviour shown in Figure 3 is reasonable, also in view of recent field findings in complex sea state conditions [45].
* = 0.63/ = 1.54 m s , the drag coefficient is levelled off and starts to decrease with a further increase in wind speed; (ii) in the considered wave age range, the mechanical drag at the sea surface decreases as the wave age increases.A possible explanation for the well-known decrease of the drag coefficient with the wind speed increase above the hurricane force (a critical threshold of about 33 m s −1 is typically assumed) is the development of a sea foam layer at the air-sea interface [21].In general, in the literature, the dependence of the drag coefficient on wave age is uncertain, because the slope of the fitted line can be influenced considerably by the narrow range of wave age [44].In addition, the presence of swell can generate complex interactions [23].Nevertheless, we believe that the behaviour shown in Figure 3 is reasonable, also in view of recent field findings in complex sea state conditions [45].
The Liu et al. [26,27] formulation-Equation ( 2)-has been implemented in WRF-Chem under the Mellor-Yamada-Nakanishi-Niino (MYNN) surface layer package and it will be described in the following Sections.

Mechanisms of Sea Spray Generation under High Wind Conditions
The estimate of sea spray emission from sea surface is a well-known modelling challenge [46].Sea spray droplets span a large range of sizes, from radii of a few nanometres to several millimetres.The past decade has recorded a remarkable progress in evaluating small (sub-micrometre) and medium-sized (radius at formation r0 < 50 μm) drop emission, both mainly due to bursting bubbles [17].On the contrary, a large scatter persists in the emission parameterizations for large drops, because the formation process is still poorly understood [17].Large drops are generated mainly when the wind shear at the surface is sufficiently high for water to be literally torn off from the surface waves, originating the so-called "spume droplets."A wind speed higher than approximately 7-11 m s −1 is typically necessary to generate spume droplets [18].Even if Fairall et al. [47], during the Spray Production and Dynamics Experiment (SPANDEX), measured spume droplets up to r0 = 700 μm, it is commonly assumed that they have radii in the range r0 = 20-500 μm with a peak centred at 100 μm The Liu et al. [26,27] formulation-Equation (2)-has been implemented in WRF-Chem under the Mellor-Yamada-Nakanishi-Niino (MYNN) surface layer package and it will be described in the following Sections.

Mechanisms of Sea Spray Generation under High Wind Conditions
The estimate of sea spray emission from sea surface is a well-known modelling challenge [46].Sea spray droplets span a large range of sizes, from radii of a few nanometres to several millimetres.The past decade has recorded a remarkable progress in evaluating small (sub-micrometre) and medium-sized (radius at formation r 0 < 50 µm) drop emission, both mainly due to bursting bubbles [17].On the contrary, a large scatter persists in the emission parameterizations for large drops, because the formation process is still poorly understood [17].Large drops are generated mainly when the wind shear at the surface is sufficiently high for water to be literally torn off from the surface waves, originating the so-called "spume droplets."A wind speed higher than approximately 7-11 m s −1 is typically necessary to generate spume droplets [18].Even if Fairall et al. [47], during the Spray Production and Dynamics Experiment (SPANDEX), measured spume droplets up to r 0 = 700 µm, it is commonly assumed that they have radii in the range r 0 = 20-500 µm with a peak centred at 100 µm [17].Although observations of breaking waves at sea suggest that spume droplets are generated at the crest of breaking waves, the spume droplets generation mechanism is much less understood than the generation of film and jet droplets, due to enormous difficulties in field experiments under high wind conditions.Troitskaya et al. [48], employing high-speed video-filming in a high-speed wind-wave flume, proposed the bag-breakup mode of fragmentation of liquid in gaseous flows to be the dominant sea spray generation mechanism in extreme winds.This regime is characterized by inflating and consequent bursting of short-lived objects, bags, comprising sail-like water films surrounded by massive liquid rims.In their experiments, both breaking projection and underwater bubble-bursting were found to be less effective.A laboratory study, the Spray Production and Dynamics Experiment (SPANDEX), has been conducted by Fairall et al. [47] at the University of New South Wales in Melbourne (Australia) in the Water Research Laboratory.The targets of SPANDEX were to investigate physical aspects of spume droplet production and transport.The obtained database shows the raw count of the droplets considering both weak and strong forcing winds on waves.
Classical schemes for the Sea Spray Source Function (SSSF) rely on the concept of whitecap coverage fraction [49], the fraction of the sea surface that is covered by white capping.In this context, SSSF depends on the wind speed solely and is valid for small and medium-size droplets.To address the dependence of the sea spray emission on the wave state and to take into account the spume droplets generation, we implemented in WRF-Chem the Wan et al. [28] parameterization, which is a revised formulation of the SSSF developed by Zhao et al. [50].Based on Wan et al. [28] the SSSF for jet, film and spume droplets is expresses by: where dF is the flux (number) of sea spray particles whose radius is in the range (r 0 , r 0 + dr), r 0 is the radius of the particle at formation, R B = u 2 * ω p υ is the wind-sea Reynolds number, u * is the friction velocity, ω p is the angular frequency spectral peak of wind waves and υ is the kinematic viscosity of air.In our case, due to the fact that during the simulation period there was no swell over the Ionian Sea, ω p was approximated by the wave peak angular frequency and using the dispersion relationship for deep water, for example, Reference [51], with the definition of the wave age parameter, β = c p /u * , we obtained [26,27]: The wind-sea Reynolds number is a non-dimensional parameter that represents the coupling effect of wind forcing and wind wave state and can be regarded as a measure of the fluid dynamical conditions at the air-sea boundary layer.Zhao et al. [50] demonstrated that the sea spray production rate for spume droplets is much better correlated with the wind-sea Reynolds number rather than with the surface wind speed alone and they concluded that spume droplets begin to be produced as R B exceeds 10 3 .

SWAN Setup
In the application discussed in the present paper, wind-wave dynamics is simulated using the SWAN model (Simulating Waves Nearshore), a third-generation wave model that computes wind-generated waves in offshore and coastal regions, coupled with WRF.All runs are performed with the COAWST package [30], following the approach recently used in Ricchi et al. [33,42].
In the two-way coupled WRF-SWAN run, WRF provides the wind field to SWAN and receives back the peak period, the significant wave height and the peak wavelength every 300 s.Starting from this information about the wave fields, WRF computes the roughness value using the Oost scheme [52].The setting for the WRF model (main parameterization schemes, large scale forcing, etc.) is the same used for the WRF-Chem model (see Section 2.5).
For energy exchange consistence and with the aim of limiting the field interpolation and reducing the noise produced at the boundaries of the atmospheric model domain, we use the same numerical grid configuration both in WRF and SWAN models.In more detail, the SWAN grid is generated starting from the WRF grid but eliminating 4 points in west-east direction and 12 points in north-south direction.As a consequence, the WRF domain consists of 250 × 200 grid points and SWAN of 246 × 188 points, with horizontal grid spacing ∆ x/y = 8 km (this is comparable with that adopted in recent studies on Medicanes [4,7,11,15,53]).WRF and SWAN runs were performed for the same time duration of 48 h.
SWAN adopts a wave action spectrum discretization with 72 equally spaced directions in order to better describe the space and time variability of the wave direction, in particular around the centre of cyclone where direction and intensity of the wind may change rapidly.The used integration time-step is 30 s.Because the presence of intense winds, like in the area surrounding the "eye" of the Medicane, is known to induce wave-breaking, the white capping option (steepness-induced wave-breaking) is activated in this simulation [54].We did not perform a wave model spin-up, a choice justified by the fact that the simulation started approximately 24 h before the event (i.e., the initial time is 0000 UTC, 25 September 2006) and during the simulation period there were no other sea state disturbances over the Ionian Sea, which actually started in calm conditions.Although a spin-up period of 24 h appeared to be sufficient in this specific case, we remind that this is a critical aspect that has to be carefully evaluated case by case.Future studies should then carefully analyse the spin-up time to be adopted, since it represents a potential source of uncertainty and error.

WRF-Chem Setup
In this work, the WRF-Chem model version 3.6.1.has been used [31] in the final stage of the coupled sequence with COAWST.Figure 2 shows the numerical domain, which covers the central Mediterranean regions with 250 × 200 grid points (∆ x/y = 8 km) and 39 vertical levels.The region of interest is shown in the internal (dotted) box.The simulation started on 25 September 2006, at 0000 UTC and ended on 27 September 2006, at 0000 UTC (48 h long).Boundary and initial conditions were provided from NCAR/NCEP Final Analysis (FNL from GFS) (ds083.2),with 1-degree horizontal resolution.
Concerning the model setting of the physical parameters, we have modified the surface layer module corresponding to the MYNN [55] scheme (sf_sfclay_physics = 5).It is based on the Monin-Obukhov similarity theory and must be used in conjunction with the MYNN Planetary Boundary-Layer parameterization scheme (bl_pbl_physics = 5).In particular, we have implemented in the MYNN surface scheme one additional option for the air-sea fluxes, which are governed by the namelist variable "isftcflx" [56].It is important to mention that the isftcflx variable only affects the exchange coefficients, by changing the values of the roughness lengths of momentum (z 0 ), sensible heat (z T ) and latent heat (z Q ).The isftcflx options for the MYNN surface layer model are depicted in Table 1.Option 5 is implemented in this work, with the purpose of connecting wind, wave state and sea spray in the roughness length calculation.The other physical schemes that are used in our simulations are the following: the RUC (Rapid Update Cycle) Land Surface Model [58] is chosen to represent the land surface interactions (sf_surface_physics = 3); the Mellor-Yamada-Nakanishi and Niino (MYNN) 2.5 level turbulent kinetic energy parameterization (bl_pbl_physics = 5) is used to describe the planetary boundary layer [55]; the short/long wave radiation effects are parameterized using the Rapid Radiative Transfer Model [59] for both short-and long-wave (ra_sw_physics = ra_lw_physics = 4); the two-moment cloud microphysics scheme of Morrison et al. [60] is used for the treatment of the microphysical processes (mp_physics = 10).

WRF/SWAN Coupling
The air-sea coupled modelling system used in this paper is based on a two-stage procedure, as shown in Figure 4.The first stage is based on the COAWST modelling system.COAWST was developed to couple the WRF model for the atmosphere side, the ROMS model for the ocean circulation and the SWAN model to simulate the wave dynamics (in the present work, only the WRF and SWAN coupling is used).In the second stage, the time varying gridded wave field, generated by the COAWST system, is supplied to the WRF-Chem model, which is able to simulate the emission and transport of sea spray aerosol through the GOCART sectional aerosol model.The other physical schemes that are used in our simulations are the following: the RUC (Rapid Update Cycle) Land Surface Model [58] is chosen to represent the land surface interactions (sf_surface_physics = 3); the Mellor-Yamada-Nakanishi and Niino (MYNN) 2.5 level turbulent kinetic energy parameterization (bl_pbl_physics = 5) is used to describe the planetary boundary layer [55]; the short/long wave radiation effects are parameterized using the Rapid Radiative Transfer Model [59] for both short-and long-wave (ra_sw_physics = ra_lw_physics = 4); the two-moment cloud microphysics scheme of Morrison et al. [60] is used for the treatment of the microphysical processes (mp_physics = 10).

WRF/SWAN Coupling
The air-sea coupled modelling system used in this paper is based on a two-stage procedure, as shown in Figure 4.The first stage is based on the COAWST modelling system.COAWST was developed to couple the WRF model for the atmosphere side, the ROMS model for the ocean circulation and the SWAN model to simulate the wave dynamics (in the present work, only the WRF and SWAN coupling is used).In the second stage, the time varying gridded wave field, generated by the COAWST system, is supplied to the WRF-Chem model, which is able to simulate the emission and transport of sea spray aerosol through the GOCART sectional aerosol model.The coupling between the two stages is ruled explicitly by the wave age parameter (β), which allows, through Equation ( 4), the calculation of the wind-sea Reynolds number (RB) and thus of the SSSF, as expressed by Equation (3a,d).
To consider the overall coupled effects of wind, waves and spume droplets we run two series of simulations utilizing the full operative sequence, that is, including COAWST and WRF-Chem.In the The coupling between the two stages is ruled explicitly by the wave age parameter (β), which allows, through Equation ( 4), the calculation of the wind-sea Reynolds number (R B ) and thus of the SSSF, as expressed by Equation (3a,d).
To consider the overall coupled effects of wind, waves and spume droplets we run two series of simulations utilizing the full operative sequence, that is, including COAWST and WRF-Chem.In the first set, the WRF-Chem model runs with the default sea surface flux parameterization (isftcflx = 0); hereinafter, this set is named CPL0.In the second set, hereinafter CPL5, the WRF-Chem model runs with the sea surface flux parameterization implemented in this work, which may be selected setting the namelist variable isftcflx = 5, as depicted in Table 1.
The outputs of the COAWST modelling system, time-averaged between 0500 and 2000 UTC on 26 September 2006, are illustrated in Figure 5a-c.This period, considering the Medicane track [6,35], corresponds to the transit of the cyclone across the northern Ionian Sea and the Adriatic Sea.The region affected by the passage of the Medicane during 26 September is identified by the internal rectangular box in Figure 3.
26 September 2006, are illustrated in Figure 5a-c.This period, considering the Medicane track [6,35], corresponds to the transit of the cyclone across the northern Ionian Sea and the Adriatic Sea.The region affected by the passage of the Medicane during 26 September is identified by the internal rectangular box in Figure 3.
The wind-sea Reynolds number (RB) is shown in Figure 5a: in the region affected by the Medicane, RB takes averaged value in the range (20-35) × 10 3 .Being RB a parameter that specifies the coupling effect between wind forcing and wind waves [50], it gives a clear indication of the fluid dynamical conditions at the air-sea boundary layer.The corresponding distribution of the average mass production rate of all droplets (r0 > 10 μm) is depicted in Figure 5b.It is clear that the highest RB values correspond to abundant production rates of such sea spray droplets, with a total mass flux in the range of (200-350) μg m −2 s −1 .Whereas this is an average value, it seems in good agreement with the estimation of Zhao et al. [50] that, analysing experimental data, found an interpolation curve for the mass production rate proportional to 10  . . Figure 5c shows the average spatial distribution of the wave age parameter,  =  / * .It is apparent that the average values in the range (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30) indicate a young sea state due to the TLC forcing.This is consistent with the rapid movement of the cyclone in the region.In conclusion, the analysis of the COAWST output reveals that in the region affected by the Medicane there is a relevant production of spume droplets, corresponding to a high value of the wind-sea Reynolds number and that a young sea state is generated by the observed strong wind.

Analysis of the WRF-Chem Output
The ECMWF analysis at 0600 UTC 26 September 2006 (Figure 6a) reveals the presence of two mean sea-level pressure minima [6,35].The first low (minimum of 1003 hPa) is situated offshore the Tyrrhenian coast, centred at approximately (40.5° N, 12.8° E).This minimum is associated with the synoptic upper level low over the Tyrrhenian Sea (Figure 1a), which is responsible for the cyclonic circulation over most of the Italian peninsula.The second minimum of 996 hPa, located over the Ionian Sea, is associated with a small scale, deep cyclone, the one which is going to get tropical features.The mean sea-level pressure field at 0600 UTC obtained in CPL0 and CPL5 runs is very similar (cfr.Figure 6b with Figure 6c).The intensity of the simulated pressure minimum is very close to the analysis, while the location of the cyclone is slightly different, since in the runs it is reproduced near the western coast of Apulia region, while in the analysis it is still near the coast of Calabria.This means that the WRF model anticipates the evolution of the cyclone by a couple of hours, as discussed The wind-sea Reynolds number (R B ) is shown in Figure 5a: in the region affected by the Medicane, R B takes averaged value in the range (20-35) × 10 3 .Being R B a parameter that specifies the coupling effect between wind forcing and wind waves [50], it gives a clear indication of the fluid dynamical conditions at the air-sea boundary layer.The corresponding distribution of the average mass production rate of all droplets (r 0 > 10 µm) is depicted in Figure 5b.It is clear that the highest R B values correspond to abundant production rates of such sea spray droplets, with a total mass flux in the range of (200-350) µg m −2 s −1 .Whereas this is an average value, it seems in good agreement with the estimation of Zhao et al. [50] that, analysing experimental data, found an interpolation curve for the mass production rate proportional to 10 −12 R 1.55 B . Figure 5c shows the average spatial distribution of the wave age parameter, β = C p /u * .It is apparent that the average values in the range (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30) indicate a young sea state due to the TLC forcing.This is consistent with the rapid movement of the cyclone in the region.
In conclusion, the analysis of the COAWST output reveals that in the region affected by the Medicane there is a relevant production of spume droplets, corresponding to a high value of the wind-sea Reynolds number and that a young sea state is generated by the observed strong wind.

Analysis of the WRF-Chem Output
The ECMWF analysis at 0600 UTC 26 September 2006 (Figure 6a) reveals the presence of two mean sea-level pressure minima [6,35].The first low (minimum of 1003 hPa) is situated offshore the Tyrrhenian coast, centred at approximately (40.5 • N, 12.8 • E).This minimum is associated with the synoptic upper level low over the Tyrrhenian Sea (Figure 1a), which is responsible for the cyclonic circulation over most of the Italian peninsula.The second minimum of 996 hPa, located over the Ionian Sea, is associated with a small scale, deep cyclone, the one which is going to get tropical features.The mean sea-level pressure field at 0600 UTC obtained in CPL0 and CPL5 runs is very similar (cfr.Figure 6b with Figure 6c).The intensity of the simulated pressure minimum is very close to the analysis, while the location of the cyclone is slightly different, since in the runs it is reproduced near the western coast of Apulia region, while in the analysis it is still near the coast of Calabria.This means that the WRF model anticipates the evolution of the cyclone by a couple of hours, as discussed in Moscatello et al. [6].Also, the distribution of the isobars is slightly different, since the isolines of 998 and 996 hPa are much closer in the simulations than in the analysis.
Atmosphere 2018, 9, x FOR PEER REVIEW 11 of 17 in Moscatello et al. [6].Also, the distribution of the isobars is slightly different, since the isolines of 998 and 996 hPa are much closer in the simulations than in the analysis.An extensive observational coverage of the event is provided in Moscatello et al. [35] and in Conte et al. [39].Using radar maps and satellite images, it was possible to identify the track of the cyclone (Figure 7, X points) from the time the small-scale cyclone crossed the Salento peninsula (the south-eastern tip of Italy) to its decline near the Gargano promontory after crossing the southern Adriatic Sea. Figure 7 shows also the simulated Medicane trajectories, that are built representing the position of the sea-level pressure minimum in the numerical domain.This is done both for CPL0 run (black-line) and CPL5 run (red-line).The Medicane tracking started at 0000 UTC, 26 September and ended at 1800 UTC of the same day.It is apparent from Figure 7 that the CPL5 run provides a better representation of the cyclone track along the Adriatic Sea.
Figure 8 shows the drag coefficient (shaded) and mean sea-level pressure (line contours) time averaged between 0500 and 2000 UTC on 26 September for CPL0 run (Figure 8a) and CPL5 run (Figure 8b).The drag coefficient at 10 m is calculated as: where Ψ is the stability function, which depends on the Monin-Obukhov length (L).The roughness formulations given by Equations ( 1) and ( 2)-used in the CPL0 and CPL5 run, respectively-are "relocated" on the drag coefficient through Equation ( 5).The average spatial pattern of the drag coefficient is different between the two simulations.This is a consequence of the wave field that, by means of the wave age parameter (β), provides an increment of the surface drag coefficient during the passage of the Medicane in the CPL5 run.This fact is particularly evident in the sub-region (42-43° N) over the Adriatic Sea, where the average wave height was 2.8 m (not shown) and <β> = 21, as already documented in Section 3.1 (Figure 5c).The sea-spray component (ω parameter) on roughness formulation expressed by Equation (2) needs greater wind speed to be dominant (>30 m s −1 ).These changes affect the sea-level pressure field, thus modify the track of the TLC as depicted in Figure 7.In particular, the pressure minimum of 998 hPa near the Gargano Promontory for the CPL5 run (Figure 8b) corresponds to a track closer to coastline and more consistent with the observations (Figure 7, red-line) in comparison with the track for the CPL0 run (Figure 7, black-line).
The common metrics on cyclone intensity, that is, the minimum sea-level pressure (SLP) and the maximum 10-m wind speed, are plotted in Figure 9a,b over the entire 48-h period of simulation.The ratio between the surface exchange coefficients of moist enthalpy (Ck) and drag (CD) is also plotted in Figure 9c.This ratio is an important scaling parameter in the theory of potential intensity [1] since the An extensive observational coverage of the event is provided in Moscatello et al. [35] and in Conte et al. [39].Using radar maps and satellite images, it was possible to identify the track of the cyclone (Figure 7, X points) from the time the small-scale cyclone crossed the Salento peninsula (the south-eastern tip of Italy) to its decline near the Gargano promontory after crossing the southern Adriatic Sea. Figure 7 shows also the simulated Medicane trajectories, that are built representing the position of the sea-level pressure minimum in the numerical domain.This is done both for CPL0 run (black-line) and CPL5 run (red-line).The Medicane tracking started at 0000 UTC, 26 September and ended at 1800 UTC of the same day.It is apparent from Figure 7 that the CPL5 run provides a better representation of the cyclone track along the Adriatic Sea.
Figure 8 shows the drag coefficient (shaded) and mean sea-level pressure (line contours) time averaged between 0500 and 2000 UTC on 26 September for CPL0 run (Figure 8a) and CPL5 run (Figure 8b).The drag coefficient at 10 m is calculated as: where Ψ m is the stability function, which depends on the Monin-Obukhov length (L).The roughness formulations given by Equations ( 1) and ( 2)-used in the CPL0 and CPL5 run, respectively-are "relocated" on the drag coefficient through Equation ( 5).The average spatial pattern of the drag coefficient is different between the two simulations.This is a consequence of the wave field that, by means of the wave age parameter (β), provides an increment of the surface drag coefficient during the passage of the Medicane in the CPL5 run.This fact is particularly evident in the sub-region (42-43 • N) over the Adriatic Sea, where the average wave height was 2.8 m (not shown) and <β> = 21, as already documented in Section 3.1 (Figure 5c).The sea-spray component (ω parameter) on roughness formulation expressed by Equation (2) needs greater wind speed to be dominant (>30 m s −1 ).These changes affect the sea-level pressure field, thus modify the track of the TLC as depicted in Figure 7.In particular, the pressure minimum of 998 hPa near the Gargano Promontory for the CPL5 run (Figure 8b) corresponds to a track closer to coastline and more consistent with the observations (Figure 7, red-line) in comparison with the track for the CPL0 run (Figure 7, black-line).
The common metrics on cyclone intensity, that is, the minimum sea-level pressure (SLP) and the maximum 10-m wind speed, are plotted in Figure 9a,b over the entire 48-h period of simulation.
The ratio between the surface exchange coefficients of moist enthalpy (C k ) and drag (C D ) is also plotted in Figure 9c.This ratio is an important scaling parameter in the theory of potential intensity [1] since the maximum wind speed and minimum sea-level pressure that a tropical cyclone may reach depend on this fraction.
Atmosphere 2018, 9, x FOR PEER REVIEW 12 of 17 maximum wind speed and minimum sea-level pressure that a tropical cyclone may reach depend on this fraction.0840 to 0940 UTC, while the recorded SLP minimum was 986 hPa at 0915 UTC.The simulated SLP is plotted in Figure 9a for the CPL0 (black line) and the CPL5 (red line) simulations.In CPL0, the cyclone crossed Salento between 0700 UTC and 0800 UTC with simulated SLP of 991 and 990 hPa, respectively.In contrast, the CPL5 run shows a better matching with the observational data, as the simulated crossing time was between 0800 and 0900 UTC and simulated SLP values of 987 and 984 hPa, respectively.The maximum 10-m wind speed is plotted in Figure 9b for the CPL0 (black-line) and the CPL5 (red-line) simulations.These two runs show a similar time evolution pattern, being the highest velocities simulated in CPL5 closer to observation.One can note that the maximum wind speed occurs in correspondence with the crossing of the cyclone over the Ionian (2-7 h) and Adriatic Sea (8-13 h).In these periods, the CPL5 run provides the ratio Ck/CD (Figure 9c), a little above the value of 0.60.This is in fair agreement with the experimental value derived from CBLAST experiment [22] that showed an average value near 0.70 for tropical storm conditions (surface wind greater than 20 m s −1 ).
Figure 10 shows the time evolution of SLP in the two simulations, taken in the grid point closest to Nardò station, in the middle of Salento peninsula, where the lowest pressure was measured, in comparison with the observations.It is apparent that both CPL0 (blue continuous line) and CPL5 (red dashed line) have a sharp drop in SLP of about 10 hPa, but, again, CPL5 reproduces a time evolution closer to measurements, both in intensity, timing and rapidity of evolution.The analysis in Moscatello et al. [35] showed that the cyclone crossed the Salento peninsula from 0840 to 0940 UTC, while the recorded SLP minimum was 986 hPa at 0915 UTC.The simulated SLP is plotted in Figure 9a for the CPL0 (black line) and the CPL5 (red line) simulations.In CPL0, the cyclone crossed Salento between 0700 UTC and 0800 UTC with simulated SLP of 991 and 990 hPa, respectively.In contrast, the CPL5 run shows a better matching with the observational data, as the simulated crossing time was between 0800 and 0900 UTC and simulated SLP values of 987 and 984 hPa, respectively.
The maximum 10-m wind speed is plotted in Figure 9b for the CPL0 (black-line) and the CPL5 (red-line) simulations.These two runs show a similar time evolution pattern, being the highest velocities simulated in CPL5 closer to observation.One can note that the maximum wind speed occurs in correspondence with the crossing of the cyclone over the Ionian (2-7 h) and Adriatic Sea (8-13 h).In these periods, the CPL5 run provides the ratio C k /C D (Figure 9c), a little above the value of 0.60.This is in fair agreement with the experimental value derived from CBLAST experiment [22] that showed an average value near 0.70 for tropical storm conditions (surface wind greater than 20 m s −1 ).
Figure 10 shows the time evolution of SLP in the two simulations, taken in the grid point closest to Nardò station, in the middle of Salento peninsula, where the lowest pressure was measured, in comparison with the observations.It is apparent that both CPL0 (blue continuous line) and CPL5 (red dashed line) have a sharp drop in SLP of about 10 hPa, but, again, CPL5 reproduces a time evolution closer to measurements, both in intensity, timing and rapidity of evolution.

Conclusions
This work represents a first attempt toward an integrated numerical system that utilises the COAWST system in combination with the chemistry package WRF-Chem.The latter is configured using the GOCART aerosol module, which explicitly considers the emission and transport of sectional sea spray aerosol particles.With this sequence, there is consistency between the wind field (peak period, significant wave height and peak wavelength) and spray induced roughness under the extreme conditions associated with hurricane winds.
The focus of this study is on the development of waves storms in wind sea conditions, due to the associated wind forcing at the sea surface under strong wind and the consequent generation of a large quantity of sea spray droplets.During strong wind conditions, large sea spray droplets are produces mainly a consequence of the mechanical tearing of the wave crests by the wind.These sea spray droplets are generally called "spume droplets" with a radius at formation in the range (20-500) μm.The presence of these spume droplets changes the basic state of the air-sea interface in the marine surface boundary layer, influencing the air-sea exchange processes of momentum and heat and modifying the sea surface drag coefficients.To describe properly these complex phenomena, the WRF-Chem model has been modified to allow the generation of spume droplets and the relative modification of the sea surface roughness length, considering the combined effects of waves and sea spray.
A numerical test is performed considering the Medicane occurred in south-eastern Italy on 26 September 2006.Two series of coupled runs have been performed: the first one (CPL0), using the default sea surface roughness formulation of Charnock [25]; the second one (CPL5), with the formulation implemented in this work, which considers both the effects of waves and spume droplets.
The CPL5 run produces a better matching to the observed cyclone characteristics, both for the track and for common cyclone metrics such as: maximum wind speed, minimum central pressure and the ratio of exchange coefficients Ck/CD.This means that the use of a coupled atmosphere-wave-

Conclusions
This work represents a first attempt toward an integrated numerical system that utilises the COAWST system in combination with the chemistry package WRF-Chem.The latter is configured using the GOCART aerosol module, which explicitly considers the emission and transport of sectional sea spray aerosol particles.With this sequence, there is consistency between the wind field (peak period, significant wave height and peak wavelength) and spray induced roughness under the extreme conditions associated with hurricane winds.
The focus of this study is on the development of waves storms in wind sea conditions, due to the associated wind forcing at the sea surface under strong wind and the consequent generation of a large quantity of sea spray droplets.During strong wind conditions, large sea spray droplets are produces mainly a consequence of the mechanical tearing of the wave crests by the wind.These sea spray droplets are generally called "spume droplets" with a radius at formation in the range (20-500) µm.The presence of these spume droplets changes the basic state of the air-sea interface in the marine surface boundary layer, influencing the air-sea exchange processes of momentum and heat and modifying the sea surface drag coefficients.To describe properly these complex phenomena, the WRF-Chem model has been modified to allow the generation of spume droplets and the relative modification of the sea surface roughness length, considering the combined effects of waves and sea spray.
A numerical test is performed considering the Medicane occurred in south-eastern Italy on 26 September 2006.Two series of coupled runs have been performed: the first one (CPL0), using the default sea surface roughness formulation of Charnock [25]; the second one (CPL5), with the formulation implemented in this work, which considers both the effects of waves and spume droplets.
The CPL5 run produces a better matching to the observed cyclone characteristics, both for the track and for common cyclone metrics such as: maximum wind speed, minimum central pressure and the ratio of exchange coefficients C k /C D .This means that the use of a coupled atmosphere-wave-aerosol model allows a non-negligible improvement in the simulation skill of the numerical system.This is obtained through the representation of heat and moisture fluxes consistent among the different modules of the numerical system and by providing a more refined, both in time and space, sea surface roughness related to wave state and sea spray aerosols.
Finally, the study is focused on the influence of wave state and sea spray on roughness length for momentum transport.However, enthalpy flux across the air-sea interface is also strongly dependent on both wave state and sea spray, as investigated for tropical cyclones and mid-latitude storms [63].The inclusion of the thermal coupling of the air-sea interaction in the modelling system is left anyway for a forthcoming study.

Figure 2 .
Figure 2. The numerical domain.The dotted box indicates the region of interest.Labels "G" and "N" denote the location of Galatina Airport and Nardò meso-network stations, respectively.

Figure 2 .
Figure 2. The numerical domain.The dotted box indicates the region of interest.Labels "G" and "N" denote the location of Galatina Airport and Nardò meso-network stations, respectively.

Figure 2 .
Figure 2. The numerical domain.The dotted box indicates the region of interest.Labels "G" and "N" denote the location of Galatina Airport and Nardò meso-network stations, respectively.

Figure 4 .
Figure 4. Schematic representation of the coupled modelling system.

Figure 4 .
Figure 4. Schematic representation of the coupled modelling system.

Figure 7 .
Figure 7. Cyclone track starting 0000 UTC 26 September 2006 of CPL0 run (continuous black line); CPL5 run (continuous red line).X points denotes the track of cyclone from radar-map images and their timing in DD/HHMM format.

Figure 7 .
Figure 7. Cyclone track starting 0000 UTC 26 September 2006 of CPL0 run (continuous black line); CPL5 run (continuous red line).X points denotes the track of cyclone from radar-map images and their timing in DD/HHMM format.

Figure 7 .
Figure 7. Cyclone track starting 0000 UTC 26 September 2006 of CPL0 run (continuous black line); CPL5 run (continuous red line).X points denotes the track of cyclone from radar-map images and their timing in DD/HHMM format.

Figure 9 .
Figure 9.Time series of the minimum sea-level pressure (a); maximum 10-m winds (b) and (c) ratio between moist enthalpy (Ck) and drag coefficient for momentum (CD).The minimum sea level pressure measured in Nardò station (986 hPa at 0915 UTC) and the maximum wind speed recorded in Galatina station (WMO station) (40 m/s at 0930 UTC) are also shown (blue asterisks) respectively in (a,b) (data from [6]) both taken along the cyclone track.

Figure 9 .
Figure 9.Time series of the minimum sea-level pressure (a); maximum 10-m winds (b) and (c) ratio between moist enthalpy (C k ) and drag coefficient for momentum (C D ).The minimum sea level pressure measured in Nardò station (986 hPa at 0915 UTC) and the maximum wind speed recorded in Galatina station (WMO station) (40 m/s at 0930 UTC) are also shown (blue asterisks) respectively in (a,b) (data from [6]) both taken along the cyclone track.

Figure 10 .
Figure 10.Time series of sea-level pressure at Nardo station for the measured data (asterisks); CPL0 run (blue continuous line); CPL5 run (red dotted line).

Figure 10 .
Figure 10.Time series of sea-level pressure at Nardo station for the measured data (asterisks); CPL0 run (blue continuous line); CPL5 run (red dotted line).

Table 1 .
The default isftcflx parameterization (0) and the additional parameterization option (isftcflx = 5) included in order to take into account wave state and sea spray effects.

Table 2 .
Modified definition of the range size bins between r min and r max as used in this work (r min and r max represent the minimum and maximum droplet radii, r eff the effective radius).

Table 2 .
Modified definition of the range size bins between rmin and rmax as used in this work (rmin and rmax represent the minimum and maximum droplet radii, reff the effective radius).