Drag Coefficient Parameterization Under Hurricane Wind Conditions

1. Introduction

Numerical wave models based on the spectral representation of the wave distribution (WAVEWATCH III [1], SWAN [2], WAM [3]) generalize the main experimental and theoretical achievements obtained over the past years of the theory of surface waves [4]. In accordance with the increasing requirements for numerical forecasts, their spatial resolution is constantly increasing, however, increasing the spatial resolution in modern models of the circulation of the atmosphere and ocean requires, in turn, improving the quality of the description of the surface layer of the atmosphere, taking into account the heterogeneity of the underlying surface, which is determined by the sea state, the presence foam or ice, complex multi-phase composition of atmospheric air in strong winds, when wave breaking leads to the formation of spray. The construction of models of such complex phenomena must inevitably be based on experimental studies. In this case, the greatest difficulty is to take into account small-scale processes near the interfaces between atmosphere and hydrosphere, which, in principle, cannot be resolved in global numerical models and are specified in parametric dependencies using the so-called bulk formulas. For the boundary layer of the atmosphere and the hydrosphere, typical examples here are the effects associated with the waves breaking, the formation of spray, foam bubbles, which is especially important in the case of high wind speeds and hurricane conditions. To date, the average intensity of emerging storms is increasing due to climate change [5]. Hurricanes are accompanied by strong winds that create large waves that are dangerous for coastal areas in terms of flooding and possible erosion of the coastline. At the same time, high-precision forecasts of the weather and the sea surface will ensure the safety of the population and minimize losses from such natural phenomena.

Taking into account the designated importance of small-scale processes in extreme weather events (storms and hurricanes), their detailed studies have recently been carried out in the framework of laboratory modeling on experimental stands. The obtained results demonstrated the great influence of marine aerosol on the Earth system, including the physics and chemistry of the atmosphere above the oceans, and marine geochemistry and biogeochemistry in general [6,7]. It affects atmospheric transparency, remote sensing and air quality. The effect of sea spray can be considered as a plausible reason for the decrease in the surface drag coefficient in hurricane-force winds. The corresponding decrease in turbulent stress is explained in [8] by the exchange of momentum between the droplets and the air flow.

In this paper, we make the estimates of the parameters of exchange processes at the boundary between the ocean and the atmosphere under hurricane conditions based on the numerical prediction of surface waves. This took into account the effect that was recently named the main mechanism for the formation of spray in strong winds, fragmentation of the "parachute" type (see [9,10,11]). The proposed parameterization was subsequently tested taking into account small-scale processes in the WAVEWATCH III (WW3) wave model.

2. Methods

WW3 is based on the numerical solution of the equation for the spectral density of wave action N in the approximation of phase averaging:

$$\frac{\partial N}{\partial t}+{\nabla }_x\dot{·\mathit{x}N}+\frac{\partial }{\partial k}\dot{k}N+\frac{\partial }{\partial \theta }\dot{\theta }N=\frac{S}{\sigma }.$$

(1)

The left hand side of the equation (1) describes the kinematics of waves, σ is the radian frequency, k is the wavenumber, θ is the wave direction. In the right hand side, there are terms that describe the wind input S_in, dissipation mainly due to wave breaking S_dis, and 4-wave nonlinear interaction of waves S_nl. The model can take into account the processes of wave attenuation due to bottom friction, three-wave interaction in shallow water, energy dissipation on bottom inhomogeneities, the influence of currents, tides, the presence of subgrid-scale islands and ice cover.

In this paper, we focus on the S_in source term and its influence on the resulting wave parameters. In general, wind wave growth S_in is described by Miles mechanism and is given as wave growth under wind input as follows [12]:

$$S_{in}=\beta N\sigma ,$$

(2)

here β is wind wave growth rate [12], it expresses wave growth by the wind input:

$$Im\sigma =\frac{1}{2}\sigma \frac{u_{*}^2}{c^2}\beta ,$$

(3)

here с is phase velocity of the wave. β coefficient depends on friction velocity u_*, which is formulated by the turbulent momentum flux:

$${\tau }_{turb}={\rho }_a〈u_x^{\text{'}}{u}_y^{\text{'}}〉={\rho }_a{u}_{*}^2,$$

(4)

here ${\rho }_a$ is air density, $u_x^{\text{'}},u_y^{\text{'}}$ - fluctuation components of wind speed.

Experimental measurement of the magnitude of the turbulent momentum flux is a difficult task. One of the most common methods is the gradient method. The gradient method uses a logarithmic law based on the Prandtl and Karman boundary layer theory for a flat plate: under conditions of neutral stratification in a layer of constant flows (where the turbulent momentum flux does not depend on height), the wind speed profile is close to logarithmic (see, for example, [13]):

$$U\left(z\right)=\frac{u_{*}}{\kappa }ln\frac{z}{z_0}.$$

(5)

Here $\kappa =0.4$ is the Von Kàrmàn constant, z₀ is the roughness parameter. By analogy with the resistance of a flat plate, the drag coefficient of the water surface is introduced, which relates the measured wind speed and the turbulent momentum flux (friction velocity):

$$C_D=\frac{{\tau }_{turb}}{{\rho }_a{U}_{10}^2}=\frac{u_{*}^2}{U_{10}^2},$$

(6)

where U₁₀ is 10 m wind speed.

The parameterization of the C_D drag coefficient is of great importance in momentum exchange in the atmospheric surface boundary layer especially at hurricane winds. Thus, to refine the wind input on waves at high (hurricane) wind speeds, it is necessary to use the parameterization of the drag coefficient, which takes into account additional effects that occur at such wind speeds, in particular, the influence of emerging spray and foam.

Wind input and dissipation in the WW3 wave model are represented by different parameterizations. The ST1 model (WAM 3) [14] is given by two empirical formulas. The first one is for the wind wave growth rate:

$$S_{in}(k,\theta )=C_{in}\frac{{\rho }_a}{{\rho }_w}max\left[0,\left(\frac{28u_{*}}{c}\cos \left(\theta -{\theta }_w\right)-1\right)\right]\sigma N\left(k,\theta \right),$$

(7)

where C_in=0.25 is a constant parameter, ${\rho }_a/{\rho }_w$ is the relation of the air and water densities, θ_w is the mean wind direction, c is the phase velocity, N(k,θ) is the is the wave action density spectrum. The second one is the parameterization of the drag coefficient proposed in [15]:

$$C_D=0.001\times (0.8+0.65U_{10}).$$

(8)

This parameterization provides a relationship between the 10 m wind speed $U_{10}$ and friction velocity $u_{*}=U_{10}\sqrt{C_D}$.

The ST6 model [16,17] also consists of the wind wave growth rate:

$$S_{in}(k,\theta )=\frac{{\rho }_a}{{\rho }_w}\sigma \gamma \left(k,\theta \right)N\left(k,\theta \right)$$

(9)

$$\gamma \left(k,\theta \right)=G\sqrt{B_n}W$$

(10)

$$G=2.8-\left(1+tanh\left(10\sqrt{B_n}W-11\right)\right)$$

(11)

$$B_n=A\left(k\right)N\left(k\right)\sigma k^3$$

(12)

$$W={\left(\frac{U_s}{c}-1\right)}^2$$

(13)

Here the omni-directional action density is obtained by integration over all directions:

$$N\left(k\right)=\int N(k,\theta )d\theta$$

(14)

The A(k) is the inverse directional spectral narrowness and is defined as:

$$A^{-1}\left(k\right)={\int }_0^{2\pi }\left[N(k,\theta )/N_{max}\left(k\right)\right]d\theta$$

(15)

$$N_{max}\left(k\right)=max\left\{N(k,\theta )\right\}.$$

(16)

for all directions $\theta \in \left[\mathrm{0,2}\pi \right]$.

$U_s$ is the scaling wind speed and, following [16]:

$$U_s=U_{10}\approx \mathsf{{\rm Y}}{u}_{*}.$$

(17)

and $\mathsf{{\rm Y}}=28$ by following [18].

The value of W takes into account the direction of the wind, considering W as the sum of the winds in the direction of the waves and the headwinds, so that they complement each other.

The drag coefficient is determined by the expression:

$$C_D={10}^{-4}\times \left(8.058+0.967U_{10}-0.016U_{10}^2\right).$$

(18)

The parameterization was proposed in [19] and takes into account the saturation and further decrease in extreme winds of water resistance at wind speeds exceeding 30 m/s. To prevent u⋆ from dropping to zero in very strong winds (U₁₀ ≥ 50.33 m/s), expression (18) was modified to give u_⋆ = 2.026 m/s.

It should be noted that the parameterization ST6 also has a negative input term to attenuate wave growth in those parts of the wave spectrum where there is a negative component of wind stress (13), in addition to positive wind input. The growth rate for headwinds is negative [20].

ST4 [21] is also the modern widely used models of wind input. This parameterization is set iteratively, based on the WAM 4 model for the "positive contribution" of the wind, but introduces corrections that reduce the drag coefficient at high wind speeds. This is achieved by reducing the contribution of wind input at high winds.

3. Results and Discussion

3.1. Simulation

In this paper, a case of Irma hurricane is considered to test different parameterizations of wind input and drag coefficient parameterizations. It took place in 30/08/2017 - 12/09/2017 in the North Atlantic Ocean. Hurricane Irma was an extremely powerful Cape Verde hurricane. Irma developed from a tropical wave near Cape Verde on August 30, then it rapidly intensified into a Category 3 hurricane on the Saffir–Simpson wind scale on August 31. The intensity fluctuated on the next days, but on September 4, Irma resumed intensifying, becoming a Category 5 hurricane on the next day. The best track of the hurricane is made by USA National Hurricane Center and is shown on the Figure 1 taken from the tropical cyclone report [22]. The red dashed line shows the area that is considered in the simulation. NDBC buoys used for the further analysis are shown on the map by yellow points.

The numerical simulations of waves under hurricane Irma conditions are performed within WW3. The modelling covered the period of 01/09/2017 – 09/09/2017. The longitudes and latitudes of the considered area varied from 14 up to 24 degrees and from 280 up to 330 degrees respectively. The frequencies ranged in the interval from 0,037 Hz up to 2,97 Hz, 24 angular directions are considered. The ETOPO1 data is used to create bathymetry and mask files. The waves developed under wind forcing from the CFSv2 reanalysis with spatial resolution of 0,205⁰ [23]. It is necessary to mention, that in this test simulation the waves developed from the calm initial conditions and there are no waves in the considered area at the beginning of the simulation. The results are obtained after few hours from the beginning of the simulation. The simulations are carried out at the IAP RAS cluster using pure MPI option of the model code.

First, default approach is considered: two state-of-the-art parameterizations of wind input (ST4 [21] and ST6 [16,17]) are tested. The resulting significant wave height distribution is obtained for two types of wind input parameterizations: ST4 and ST6. The model reproduces the rotating and motion of the hurricane, and the asymmetry of the wave field distribution can be seen. On the Figure 2 the typical phases of a hurricane development and the increase of its intensification are shown. ST6 parameterization gives values higher than ST4 on approximately 15%. On September 6, the following hurricane Jose is also obtained.

3.2. Estimation of the air-sea exchange processes parameters

On the basis of the obtained in the simulations data, the parameters of exchange processes at the ocean - atmosphere boundary layer under hurricane conditions are estimated. Using the results of the wind and wave parameters obtained in the calculations, we estimated the surface drag coefficient and the effects associated with the formation of sea spray. The small-scale effect of "parachute" type fragmentation is taken into account [9,10,11]. The mechanism is the inflation and subsequent rupture of short-lived, sail-like elevations of the water surface. This process is similar to the “parachute” type fragmentation of liquid drops and jets in gas flows.

Based on the obtained data, a spray generation function is proposed for the parachute-type rupture mechanism, which depends on the wave fetch and satisfies both laboratory and field conditions, using the dimensionless Reynolds number (Re) [24]:

$${\mathit{R}\mathit{e}}_{\mathit{B}}=\frac{{{\mathit{u}}_{\mathrm{*}}}^2}{{\mathit{\sigma }}_{\mathit{p}}\mathit{\nu }},$$

(19)

where u_* is wind friction velocity, σ_p is peak frequency of surface waves, ν is the air kinematic viscosity. The function that best describes the experimental data is selected, it shows that the exact number of "parachutes" can be determined by the following equation (see Figure 3 a):

$$〈{\mathrm{N}}_{\mathrm{p}}〉={\mathrm{M}}_0{{\mathit{R}\mathit{e}}_{\mathit{B}}}^{2/3}exp\left(-\frac{{\mathrm{M}}_1}{{{\mathit{R}\mathit{e}}_{\mathit{B}}}^{2/3}}\right),$$

(20)

where M₀=2.08, M₁=972.

Taking into account the statistics of events associated with the formation of spray, the contribution of spray to the processes of exchange between the atmosphere and the ocean is estimated (see [11,25]). In particular, this makes it possible to explain the unusual nonmonotonic dependence of the sea surface drag coefficient on wind speed (see Figure 3b).

The curve in Figure 3b is obtained for a certain wave age parameter, defined as the ratio of the 10 m wind speed to the phase velocity of the waves corresponding to the peak in the frequency spectrum of surface waves c:

$$\mathit{\Omega }=\frac{{\mathit{U}}_{10}}{\mathit{c}},$$

(21)

The surface drag coefficient is also estimated for the wave parameters calculated in the WW3 model for Hurricane Irma. On Figure 4 the examples of the distribution of parameters that determine the formation of spray and other exchange parameters associated with it are shown. Figure 4a shows the wind speed and direction distribution in Hurricane Irma according to the CFSv2 reanalysis. On Figure 4b, the frequency of the spectrum peak distribution calculated in the WW3 model is shown. On Figure 4c, the distribution of the surface drag coefficient estimated for the parameters in Figure 4 a, b, taking into account the contribution of sea spray and foam is shown. A decrease in the surface drag coefficient in the area of strong winds is demonstrated. The pronounced asymmetry of the distribution of surface drag coefficient in the direction of the storm center motion reflects the asymmetry of the surface wind and the wave distribution.

Thus, the estimates allow us to conclude that it is important to take into account the effect of spray and foam in the numerical simulation wave parameters. To take into account the effect of spray and foam in the parameterization of drag coefficient and, subsequently, in the wave model, it is convenient to obtain a formula expression for CD. Such an expression is obtained [30] in the framework of numerical experiments with the Lagrangian stochastic model (LSM), the drag coefficient depends on 10 m wind speed and wave age parameter and it is different for two ranges of wind speeds:

$$\left\{\begin{array}{c}15m/s<{\mathbf{U}}_{10}<30m/s\\ {\mathrm{C}}_{\mathrm{D}}=\left(0.00231-0.000308\mathrm{\Omega }+0.00011{\mathrm{\Omega }}^2-1.29527\mathrm{E}-005{\mathrm{\Omega }}^3\right)+\\ (-4.481E-005+4.68677E-005\mathrm{\Omega }-1.6789E-005{\mathrm{\Omega }}^2+\\ 1.9389E-006{\mathrm{\Omega }}^3){\mathrm{U}}_{10}+\\ \left(1.31\mathrm{E}-006-1.678\mathrm{E}-006\mathrm{\Omega }+6.12\mathrm{E}-007{\mathrm{\Omega }}^2-6.8\mathrm{E}-008{\mathrm{\Omega }}^3\right){{\mathrm{U}}_{10}}^2+\\ (1.48E-008+1.09E-010\mathrm{\Omega }-7.922E-010{\mathrm{\Omega }}^2){{\mathrm{U}}_{10}}^3\end{array}\right.,$$

(22)

$$\left\{\begin{array}{c}30m/s<{\mathbf{U}}_{10}<70m/s\\ {\mathrm{C}}_{\mathrm{D}}=\left(-0.00716-0.00435\mathrm{\Omega }+0.00065{\mathrm{\Omega }}^2\right)+\\ +(0.0005849+0.000294\mathrm{\Omega }+4.55E-005{\mathrm{\Omega }}^2){\mathrm{U}}_{10}+\\ \left(-1.073\mathrm{E}-005+5.82\mathrm{E}-006\mathrm{\Omega }+9.32062858\mathrm{E}-007{\mathrm{\Omega }}^2\right){{\mathrm{U}}_{10}}^2+\\ (6.45E-008-3.79E-008\mathrm{\Omega }-6.316E-009{\mathrm{\Omega }}^2){{\mathrm{U}}_{10}}^3\end{array}\right.$$

(23)

The feature of the obtained expression for C_D includes the dependence not only on the wind speed, but also on the wave age, which is especially important at high wind speeds. Wave age varies in time and the values are completely different for the short waves that are in the beginning of the simulation (Ω=1 ÷ 1.5) to hurricane waves (Ω is up to 3) and then for the swell (Ω<1). On Figure 5 the change in the value of the wave age at one of the points in space where the hurricane passed is displayed. It can be seen that at first the wave age corresponds to the presence of short waves, however, from the moment the hurricane wind approaches, a sharp change in the peak frequency of the wave occurs, since in the two-peak model for describing the wave spectrum long waves begin to exceed short waves.

For simulation in the WW3 model, the parameterization (22, 23) is used together with the parameterization of the wave growth rate β with respect to ST1 and ST6. For the initial range of wind speeds (up to 18.67 m/s and 17.36 m/s, respectively), the drag coefficient parameterization (8) and (18) respectively is used.

It should be noted that the default ST1 parameterization with drag coefficient parameterization (8) cannot be used in hurricane conditions, precisely because the C_D parameterization is not applicable to such conditions, does not have a decrease or stabilization of the C_D value. If the proposed parameterization (22, 23) is used, this problem is removed, and the use of parameterization is permissible. On Figure 6 a, b a comparison of the parameterizations (8) and (18) and the new parameterization (22, 23) correspondingly is given. The new parameterization in Figure 6 is constructed with the value of wave age Ω=3, which is typical for hurricane conditions. The new expression for C_D has a maximum shifted to higher wind speeds than, for example, (18). It corresponds to [31], where it is shown that the C_D parameterizations obtained in the open ocean have an inflection at much higher wind speeds than the parameterizations obtained in the coastal zone. The new parameterization, shown in Figure 6b, has a larger area under the curve than C_D ST6 (18), in the case of high wind speeds. However, it is important to note that a change in the wave age will change the C_D values; as the wave age decreases, it will increase. However, in the range of wind speeds 20 m/s<U₁₀<30 m/s, expression (18) provides larger values of C_D than (22, 23). Nevertheless, under hurricane conditions significant wave heights calculated using the new parameterization can be expected to be higher than those calculated using the old ones.

3.3. Implementation of new drag coefficient parameterization

To use the CD parameterization in the form as in Figure 6, the program codes w3flx1md and w3flx4md are modified for the ST1 and ST6 parameterizations, respectively. In the w3flx1md module, in the speed range U10<18.67 m/s the expression (8) is used, in the w3flx4md module, in the speed range U10<17.36 m/s, the expression (18) is used. Further, for both parameterizations, the expression (22, 23) was specified, and (22) was performed for ST1 in the speed range 18.67m/s<U10<30 m/s, for ST6 - in the speed range 17.36 m/s<U10<30 m/s (see Figure 6). Using these modified blocks, simulations are performed for the Hurricane Irma test case. The model performed well, it also reproduces the rotating and motion of the hurricane and the asymmetry of the wave distribution. On Figure 7 simulations using the ST1 (left) and ST6 (right) parameterizations are shown, each of which is applied with the CD spray specified by the formula (22, 23), respectively. This comparison allows us to see the difference in the way β is set. Waves simulated under the hurricane using ST6, where the parameter β is set in a more complex way than in ST1, taking into account, for example, a negative wave growth rate in a headwind, are better grouped and have a typical structure. For the ST1 case, some "blurring" of surface waves away from the main waves under the action of a hurricane wind is observed. It does not correspond Figure 2, where the simulation using the default models ST4 and ST6 are shown.

The simulated significant wave height of surface waves induced by the hurricane is compared with the NDBC buoys data. We examine significant wave height (Hs) and 10 m wind speed at three NDBC buoys. The buoy 41040 is located at (307.28, 14.45); 41043 is located at (295.17,21.12); 41044 is located at (301.26, 21.71). Buoys location is shown in the Figure 1. The NDBC buoys have a frequency range from 0.02 to 0.485 Hz. The wave instrumentation and data acquisition of the NDBC buoys are described in [32].

First, the wind measured on the buoys and the wind in the reanalysis that is used in the wave simulation at the same points are compared. The CFSv2 reanalysis reproduces well the shape of the wind recorded by the buoy but it overestimates the wind speed measured on the buoys especially at high wind speeds. For the buoy 41044, it is a slight overestimation, for the buoy 41040 it's more noticeable because it has rather short peaks at its time series. (Figure 8) The most striking difference is shown for the buoy 41043, it records the highest wind speeds and it is overestimated on approximately 35% (Figure 7 b).

The comparison of the simulated significant wave height of surface waves and the observational ones is performed. In the beginning of the simulation, the calculated values underestimate the buoys measurements due to the simulation assumption, according to which the simulation is carried out from the calm initial conditions, respectively, not all the swell coming from the outside is considered. Then, when the hurricane winds are coming, the values of Hs in the simulations increase rapidly whereas the buoy data shows smooth growth. The predicted extreme values under high speed wind conditions have close results for all the used parameterizations. In the Figure 8c the good correlation between the measured and observed significant wave heights are obtained. The measured wind on the buoy 41044 (Figure 8c) matches the CFSv2 reanalysis wind, and the simulated values of Hs using ST4 and ST6 default parameterizations fit the measured one but it doesn't describe the sharp increase. At the same time, ST1 and ST6 with the CD spray fit the sharp increase but overestimate the decrease of Hs after the sharp increase. ST6 with the CD spray describes the result more variably in time, which corresponds to the measurements, in contrast to the smooth curve obtained as a result of the ST1 with the CD spray simulation. In the Figure 8a the described tendency persists, however, taking into account the difference in the measured wind speed and the reanalysis allows us to conclude that the wind correction will lead to a better fit of ST6 with the CD spray to the measurements. Besides, ST6 with the CD spray shows good temporal variability consistent with measurements. In the Figure 8b the biggest difference in the measured wind speed and the reanalysis is observed. It leads to the overestimation of the measured significant wave height by all the parameterizations. ST6 with the CD spray also shows good temporal variability for this buoy. This allows us to conclude about the importance of small-scale processes when modeling waves at hurricane winds and the prospects for their inclusion in modern numerical wave models.

4. Conclusions

The estimates of the parameters of exchange processes at the boundary of the ocean and the atmosphere under hurricane wind conditions are performed. For the estimates, the numerical simulation within WAVEWATCH III wave model is performed. First, the wind and wave parameters obtained in the simulation are used to estimate the drag coefficient parameterization including the exact measured number of "parachutes" in the spray [9,10,11]. Then the parameterizations of wind input are tested within WAVEWATCH III wave model. They include default ST4 and ST6 parameterizations and the ST1 and ST6 parameterizations used together with the proposed drag coefficient parameterization [30]. This parameterization is implemented in WAVEWATCH III wave model on the basis of w3flx1md and w3flx4md modules to use in ST1 and ST6 parameterizations respectively. The obtained results are compared with the NDBC buoys data. ST6 with the C_D spray shows good temporal variability for all the considered buoys. Also, the proposed parameterization demonstrated good tracking of sharp peaks of significant wave height. It can be explained by the fact that the parameterization depends not only on the wind speed, but also on the age of the waves. This is especially important in hurricane conditions. The variability of wave age during the calculation in the presence of a hurricane is demonstrated: short waves in the beginning of the simulation have wave age values in the range Ω =1 ÷ 1.5, and then for the waves under the hurricane wind speed it changes to the values up to 3. Thus, the importance of taking into account small-scale processes at the ocean and the atmosphere boundary layer, such as spray and foam, in numerical models of wind waves is demonstrated. Further increase in spatial resolution and accuracy of the simulation shows the prospect of such studies and the need to refine the parameterizations for the case of extreme weather conditions.