High-Wind Drag Coefficient Based on the Tropical Cyclone Simulated With the WRF-LES Framework

The numerical simulation of tropical cyclones has been increasingly conducted using the advanced Weather Research and Forecast (WRF) model with the large-eddy simulation (LES) technique. Given the importance of the boundary wind profile for the vertical exchange of horizontal momentum between the atmosphere and the ocean, the drag coefficient was evaluated in the numerical simulation with the WRF-LES framework at the finest horizontal grid spacing of 37 m. In the absence of the TC–ocean interaction, the drag coefficient derived from the simulated wind profile does not show the leveling off or decrease in the strong wind conditions. The drag coefficient increases with the increasing near-surface wind speed and agrees well with the extrapolation of the Large and Pond formula in the strong wind conditions. It is suggested that the boundary wind structure simulated with the LES technique may be unrealistic when the TC–ocean interaction is not fully considered.


INTRODUCTION
Numerical simulation has been a powerful tool for studying tropical cyclones (TCs). In numerical models for TC simulation, the drag coefficient (C d ) is important for determining the vertical exchange of horizontal momentum between the atmosphere and the ocean. Great effort has been made to evaluating the drag coefficient through observational analysis (Powell et al., 2003;Jarosz et al., 2007;Vickery et al., 2009;Bell et al., 2012) and theoretical studies (Donelan et al., 2004;Donelan, 2018), especially as wind speed is stronger than 30 ms −1 (French et al., 2007;Jarosz et al., 2007;Bell et al., 2012;Edson et al., 2013). It is generally believed that the drag coefficient in high-wind condition does not increase linearly with surface wind speed although the mechanisms responsible for the reduction of the drag coefficient have not been fully understood (Donelan, 2018). However, relatively few studies have been conducted to examine the drag coefficient in the numerical simulation, while the fine-scale (less than 1,000 m) features have been explicitly simulated over the past decade (Zhu, 2008;Rotunno et al., 2009;Zhu, 2013;Green and Zhang, 2015;Wu et al., 2018, Wu et al., 2019Jiang et al., 2020;Zheng et al., 2020;Zhou et al., 2020).
The large-eddy simulation (LES) technique has been incorporated into the Weather Research and Forecast (WRF) model (Mirocha et al., 2010). The WRF-LES framework has been used to conduct TC simulation with a horizontal grid spacing less than 1,000 m. For example, (Zhu, 2008) simulated the fine-scale structures in the TC boundary layer with the 300-m and 100-m spacings. Rotunno et al. (2009) conducted idealized experiments on the f-plane and found a sharp increase in randomly distributed fine-scale turbulent eddies or gusts when the horizontal grid spacing was decreased from 185 to 62 m. Recently, the WRF-LES framework has been used to simulate fine-scale structures such as horizontal rolls and tornado-scale vortex in the TC boundary layer, and extreme wind gusts have been simulated (Ito et al., 2017;Wu et al., 2018;Wu et al., 2019;Cécé et al., 2021).
With the LES technique, the energy-producing scales of threedimensional atmospheric turbulence are explicitly resolved, while the smaller-scale portion of the turbulence spectrum is removed from the flow field using a spatial filter. In other words, the vertical exchange of horizontal momentum between the atmosphere and the ocean, which is otherwise parameterized without the LES technique, is partially resolved in the WRF-LES framework. Many studies have demonstrated that the WRF-LES framework can reasonably well simulate the boundary layer structure of TCs; it is unknown whether the drag coefficient in simulated TCs is reasonably comparable to the observation (e.g., Powell et al., 2003). These fine-scale structures are closely associated with the wind profile and the vertical momentum transfer in the boundary layer (e.g., Zhu, 2008, Zhu, 2013Ito et al., 2017;Wu et al., 2018;Wu et al., 2019); it is necessary to examine the drag coefficient based on the TC simulation with the WRF-LES framework.

THE WIND PROFILE METHOD
Following Powell et al. (2003), we use the wind profile method to estimate the drag coefficient, which is based on the log-profile of wind at the bottom of the marine atmospheric boundary layer. The wind profile in neutral conditions can be written as where U is the wind speed at altitude z and the Karman constant k is set to be 0.41. The remaining variables in Eq 2, namely, the friction velocity U p and the roughness length z 0 can be estimated through logarithmic linear regression. Then C d can be derived from where U 10 and ρ are the wind speed and air density at the 10-m altitude.

DATA
For comparison with the simulation, we use wind profiles from the GPS dropsondes deployed in 120 TCs over 17 years (Wang et al., 2015). There are over 12,000 quality-controlled data profiles. In this study, we used 1,003 profiles that were measured in high wind. Examination indicates that most of the GPS real dropsondes (RDs) were deployed at the altitude of ∼6,000 m or below within a radius of less than 30 km from the TC center. The individual wind profiles are categorized into five groups based on the mean boundary layer wind (MBW), which is defined as the mean wind speed below 1,200 m. The five highwind groups correspond to MBW in the ranges of 20-29 ms −1 (531 profiles), 30-39 ms −1 (247 profiles), 40-49 ms −1 (104 profiles), 50-59 ms −1 (68 profiles), and 60-69 ms −1 (53 profiles).
The simulation data used in the study are from a numerical simulation conducted with the WRF-LES system (Wu et al., 2018;Wu et al., 2019). The simulated TC evolves over the open ocean in the large-scale background of Typhoon Matsa (2005). Note that there was no TC-ocean interaction in this simulation. Six twoway interactive domains are embedded in the outermost domain. The finest grid spacing is 1/27 km (or about 37 m) and the model consists of 75 vertical levels with 19 levels below 2 km. The domains with the grid spacing less than 1 km move with the TC. In this study, a 22-min subset at 3-s intervals from the 30th hour of wind field data (in steady state) is used (Wu et al., 2018). The data cover an area of 90 × 90 km 2 in the inner core region (the eye and eyewall). Figure 1 shows the radius-height cross section of the azimuthal mean tangential and radial winds of the simulated TC. The tangential and radial winds are averaged over 26-36 h with the data in the 1/9-km domain. The maximum tangential wind occurs around the height of 400 m and the radius of maximum wind is ∼27 km. Comparing to Zhang et al. (2011) and Ren et al. (2020), the altitude of the boundary layer jet in the simulated TCs is lower than the observation. The strongest boundary layer inflow of more than 16 ms −1 can be found at ∼32 km near the surface, and the outflow above the inflow is much weaker. In the composite of Zhang et al. (2011) and Ren et al. (2020), the strongest inflow is between 1 and 2 times RMW. Ren et al. (2020) found that the strongest inflow is 20.7 ms −1 for major hurricanes. We can see that the wind structures of the wind in the TC boundary layer are well comparable to the observation. Two types of wind profiles are derived from the 22-min wind data at 3-s intervals. The first kind of profiles is simultaneously collected in a vertical column of the wind field. To mimic the real GPS dropsonde profiles and make more direct comparison, we construct the second kind of dropsonde profiles by simulating the dropsonde falling in the simulation. The simulated dropsonde (SD) falls only in response to gravity and drag force from the wind, and the latter is given by Hane (1975), Hock and Franklin (1999), and Stern and Bryan (2018) where V r → denotes the relative speed of dropsonde to the wind, A p is the area of the parachute, ρ a is the density of surrounding air, and C ds is the drag coefficient of the dropsonde. Stern and Bryan (2018) noted that the horizontal velocity of the dropsonde is almost the same as the surrounding wind speed. Therefore, the wind drag is almost only in the vertical direction. Like the GPS dropsondes in the real TCs, the simulated dropsonde is released at the altitude of 6,000 m during the first 10 min of the simulated data. The horizontal locations and time are randomly selected in the eyewall (17.5-30 km). The drag force in Eq 4 is set to be the same as the weight of the dropsonde, leading to the dropsonde falling at the speed of the real GPS dropsonde. Then the dropsonde locations at various vertical levels can be obtained by integrating a trajectory model. A total of 3,895 dropsondes are effectively constructed for high-wind conditions. Similar to the GPS dropsondes, the individual wind profiles are also categorized into five categories. The five high-wind groups correspond to MBW in the ranges of 20-29 ms −1 (316 profiles), 30-39 ms −1 (565 profiles), 40-49 ms −1 (660 profiles), 50-59 ms −1 (1,064 profiles), and 60-69 ms −1 (1,290 profiles). The second type of wind profile is simultaneously collected in the vertical. Unlike the RD and SD profiles, the wind data at various altitudes are simultaneous and there are no horizontal drifts. The horizontal locations and time are also selected randomly. There are 7,151 wind profiles in high-wind conditions in the ranges of 20-29 ms −1 (495 profiles), 30-39 ms −1 (1,138 profiles), 40-49 ms −1 (1,335 profiles), 50-59 ms −1 (2,448 profiles), and 60-69 ms −1 (1735 profiles). For convenience, we call the second type the idealized dropsonde (ID) profiles. The purpose of the ID profiles is to evaluate the sampling error in the SD wind profiles.

THE ESTIMATED DRAG COEFFICIENT
In this study, we linearly interpolate wind profile (horizontal component) with height onto every level that is 20 m apart and starts with 10 m. For the model output, the height is calculated from the geopotential height. For the observational data, the GPS height is used when estimating the altitude of the dropsonde. For simplicity, following Powell et al., 2003, we set the boundary layer to be 1200 m thick. Taking the categories of 20-29 ms −1 , 40-49 ms −1 , and 60-69 ms −1 as an example, Figure 2 shows the mean horizontal wind and their variance of ID, SD, and RD from the same MBW bins. The fitting of the logarithmic linear regression is also plotted in Figure 2.
Since the dropsondes data are known to be heavily contaminated at the 10-m level (Jarosz et al., 2007), they were not included in the regression. Figure 2 indicates that the simulated wind profiles (SD and ID) generally follow the log-profile, as shown in Figure 2A. The boundary layer jet in the model as shown in Figures 2B,C is lower than that observed, which is likely a result of the simulated TC being stronger than most TCs observed using dropsondes (Kepert and Wang, 2001). The drag coefficient derived from the three types of wind profiles is shown in Figure 3. We used a bootstrap method to calculate the error bar. That is, a sampling with replacement of the same size as the dataset used was conducted. We calculated drag coefficients from each sample and used their 2.5 and 97.5% percentile as the boundaries of the 95% confidence interval. The drag coefficient based on the RD profile is between 1.0 × 10 −3 and 2.0 × 10 −3 and peaks for the category of 40-49 ms −1 . For comparison with the previous studies, note that horizontal coordinate is the wind speed at 10 m. The resulting drag coefficient agrees very well with the result in Powell et al. (2003).
For the SD and ID profiles, as shown in Figure 3, the derived drag coefficient increases with the increasing MBW speed, although it is close to the drag coefficient of the RD for the category of 20-29 ms −1 . One may be wondering about the large spread of the derived drag coefficient from the RD profiles. One possible reason is that the RD profiles were from different TCs, while the SD and ID profiles are from the same TC. An important feature of the derived drag coefficient from the SD and ID profiles is that it increases linearly with the MBW or near-surface wind speed. For the category of 60-69 ms −1 , it can be as large as 4 × 10 − 3 , almost twice as much as the typical observation. It is interesting that the drag coefficient in the simulation is in good agreement with the extrapolated classic Large and Pond formula (Large and Pond, 1981). Note that the formula of Large and Pond (1981) is only valid for wind speed less than 25 ms −1 .
As mentioned above, the wind data of the SD profiles at various altitudes are not simultaneous and there are horizontal drifts. In addition, the sampling bias may be due to the boundary layer inflow or the presence of small-scale coherent structures, such as tornadoscale vortices and boundary layer rolls since the dropsondes can be repelled from strong updrafts in the boundary layer. For most of the MBW bins as shown in Figure 3, the differences in the drag coefficient between the SD and ID profiles can be 8% in the high-wind condition. Given the same MBW, the drag coefficient of the SD profiles has a smaller u 10 and larger C d . It is suggested that the sampling bias has some influence on the derived drag coefficient, but the sampling bias has little influence on the increasing drag coefficient with the increase of the MBW speed.

SUMMARY
In this study, based on the high-resolution simulation of the TC, the drag coefficient is calculated with the wind profile method (Powell et al., 2003). The numerical experiment was conducted over the open ocean using the WRF-LES model at the finest grid spacing of 37 m. While the simulated drag coefficient is similar to the observation in the category of 20-29 ms −1 , it does not show the leveling off or decrease in the strong wind conditions, likely due to the fact that in LES framework cannot take the complicated response from the sea surface, including sea foam, wave-breaking, and sheltering effect. For the strong wind conditions, the drag coefficient can be as large as 4 × 10 − 3 , almost twice as much as the typical observation. It is interesting that the drag coefficient in the simulation is in strong agreement with the extrapolated classic Large and Pond formula (Large and Pond, 1981). Since it has been demonstrated that the drag coefficient in high-wind condition does not increase linearly with surface wind speed (Powell et al., 2003), this study suggests that the boundary wind structure of the TC simulated with the LES technique may be unrealistic when the TC-ocean interaction is not fully considered.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
LW conducted analysis and writing. WJ conducted analysis and writing. QL conducted the numerical experiment.