Fluxes of momentum and moist enthalpy across the air–sea interface are believed to be one of the most important factors in determining tropical cyclone intensity. Because these surface fluxes cannot be directly resolved by numerical weather prediction models, their impacts on tropical cyclones must be accounted for through subgrid-scale parameterizations. There are several air–sea surface flux parameterization schemes available in the Weather Research and Forecasting (WRF) Model; these schemes differ from one another in their formulations of the wind speed–dependent exchange coefficients of momentum, sensible heat, and moisture (latent heat). The effects of surface fluxes on the intensity and structure of tropical cyclones are examined through convection-permitting WRF simulations of Hurricane Katrina (2005).
It is found that the intensity (and, to a lesser extent, structure) of the simulated storms is sensitive to the choice of surface flux parameterization scheme. In agreement with recent studies, the drag coefficient CD is found to affect the pressure–wind relationship (between minimum sea level pressure and maximum 10-m wind speed) and to change the radius of maximum near-surface winds of the tropical cyclone. Fluxes of sensible and latent heat (i.e., moist enthalpy) affect intensity but do not significantly change the pressure–wind relationship. Additionally, when low-level winds are strong, the contribution of dissipative heating to calculations of sensible heat flux is not negligible. Expanding the sensitivity tests to several dozen cases from the 2008 to 2011 Atlantic hurricane seasons demonstrates the robustness of these findings.
Tropical cyclones (TCs) can pose great threats to lives and property. To minimize such losses, all interested parties (governments, businesses, and individual residents) need accurate forecasts of TC position and intensity, preferably several days before the onset of hazardous weather. Although our understanding of TCs has improved considerably over the past several decades in terms of theory (e.g., Emanuel 1986, 1995b), observations (Marks et al. 2008), and modeling (e.g., Bryan 2012), progress in TC prediction has been uneven. Since 1990, the skill of position (track) forecasts has increased significantly, whereas the accuracy of intensity forecasts has not changed (Rappaport et al. 2009). The reason for this uneven progress is straightforward: TC track is mainly determined by large-scale steering flows (tropical easterlies and midlatitude westerlies) that are increasingly better resolved by global numerical weather prediction (NWP) models, whereas TC intensity is also influenced by inner-core dynamics and smaller-scale processes (such as fluxes across the air–sea interface, moist convection, radiative transfer, and cloud microphysics) that are often poorly resolved or parameterized by global and even regional models (e.g., Zhang et al. 2011). Additionally, the processes that govern intensity are inherently more chaotic and less predictable than those that govern track (e.g., Zhang and Sippel 2009).
Research groups around the world are hoping to achieve a fundamental breakthrough in the ability to forecast TC intensity, which would greatly benefit society. Many of these groups use the regional Weather Research and Forecasting (WRF) Model (Skamarock et al. 2008), which allows users to choose from an assortment of physics parameterization schemes. These different options are quite valuable because they enable researchers to isolate the effects of a certain physical process on a phenomenon of interest. Unfortunately, the freedom of choice also comes at a cost: some physics options (if not outdated altogether) may be inappropriate for the phenomenon of interest. As an example, Kepert (2012) found that TC simulations were quite sensitive to the choice of planetary boundary layer scheme, with one class of parameterization—the “Bulk” or “Hi-Res” found in the fifth-generation Pennsylvania State University–National Center for Atmospheric Research (PSU–NCAR) Mesoscale Model (MM5, the precursor to WRF)—so flawed that it “should not be used.” Furthermore, users must remember that all schemes are limited by the inherent uncertainties in the parameterizations.
In the present work, we investigate parameterizations of momentum and moist enthalpy [which includes both sensible and latent heat; see Eq. (1) in Emanuel (1995b)] fluxes across the air–sea interface, and their effects on the simulation of TCs. As discussed below, it is widely believed that a key determinant of simulated TC intensity is the bulk exchange coefficients used for surface flux calculations; the values of these coefficients are highly uncertain in strong winds over water. There have been multiple investigations of the effects of wind speed–dependent exchange coefficients on TCs (Braun and Tao 2000; Bao et al. 2000, 2002, 2012; Moon et al. 2007; Davis et al. 2008; Nolan et al. 2009a,b); of these studies, only Bao et al. (2002) focused solely on exchange coefficients, and, in that case, only sensible heat was considered. To our knowledge, this study is the first to undertake a thorough examination of different (wind speed dependent) formulations (controlled by the “isftcflx” option in the WRF namelist file) for both momentum and moist enthalpy exchange coefficients in numerically simulated TCs. Our goals are twofold:
Document the differences between distinct classes of isftcflx options in the most recent release of the Advanced Research core of WRF (WRF-ARW, version 3.4.1), including the historical background and rationale wherever possible.
Investigate how different flux options affect the structural characteristics of a simulated TC (beyond the common metrics of maximum wind speed, minimum central pressure, and track), and provide physical explanations of these differences.
The remainder of the paper is organized as follows. Section 2 provides a background discussion of air–sea momentum and moist enthalpy fluxes, as well as their impact on TC intensity and structure. Section 3 details the experimental setup. Sections 4 and 5 present results in terms of temporal evolution and radius–height structure, respectively. A general discussion can be found in section 6, followed by concluding remarks in section 7.
Tropical cyclones develop and strengthen over warm ocean waters (e.g., Gray 1968) and dissipate over land or cold water. Obviously, large fluxes of sensible heat and water vapor from the surface are a necessary (but insufficient) condition for TC genesis, intensification, and maintenance. In a very general sense, TCs may be thought of as Carnot engines, with the energy source being the aforementioned surface moist enthalpy fluxes and the sinks being dissipation in the boundary layer and in the upper-level outflow at large radii (Emanuel 1986, 1988). In this section, we review a theoretical framework that relates TC intensity to surface fluxes, discuss the processes that govern these fluxes over the ocean, and detail how these processes are parameterized in WRF-ARW.
a. Theoretical support for the importance of surface fluxes on TC intensity
The well-known potential intensity (PI) theory of Emanuel (1986)—which assumes a steady-state TC in gradient and hydrostatic balance—gives analytic solutions for the maximum tangential wind speed Vmax and minimum sea level pressure (SLP) Pmin. Emanuel (1995b) derived the nondimensional formulas [normalized by a set of scaling parameters, see Table 1 of Emanuel (1995a)] for Vmax and Pmin:
where Ck and CD are the surface exchange coefficients of moist enthalpy and drag (momentum), respectively (see the appendix); r0 is the (normalized) outer radius of the TC at which the surface wind vanishes; H is the ambient relative humidity; Ts and To are the surface and outflow temperatures, respectively; χs is the background entropy deficit (with respect to the ocean) of the subcloud layer; and Rd is the gas constant for dry air. From (1) and (2), the ratio Ck/CD clearly is an important factor in the pressure and wind of a mature storm. Indeed, some of the earliest numerical simulations of TCs (Ooyama 1969; Rosenthal 1971) found intensity to be very sensitive to changes in CH and CQ (the exchange coefficients for sensible heat and latent heat, respectively), as well as CD; this sensitivity has also been found in more recent modeling studies (e.g., Braun and Tao 2000; Bao et al. 2002; Bryan 2012).
b. Physics of exchange coefficients
The exchange coefficients are a function of several variables, including wind speed, underlying surface characteristics (ocean and different types of land), roughness elements, and atmospheric stability. Because our focus is on TCs, we only consider the bulk exchange coefficients of the sea surface.
1) Drag coefficient
Although the sea surface does not have any of the static features found on land (such as terrain), there are transient waves that temporarily change the roughness characteristics. Above wind speeds of 5 m s−1, surface gravity waves are the main roughness feature (Fairall et al. 2003) and, on their own, cause CD to increase with wind speed (Smith 1988). Because there were no direct flux measurements over the open ocean for wind speeds above 22 m s−1 until the early 2000s (Black et al. 2007), researchers were forced to extrapolate data from low wind speeds.
Once measurements near the sea surface in high wind speed conditions began to be taken, it became clear that extrapolation for the drag coefficient was inaccurate. Using GPS dropwindsonde data, Powell et al. (2003) reported that CD began to level off and then decrease for wind speeds above 33 m s−1. Very similar results were found in a laboratory tank experiment by Donelan et al. (2004); those authors noted that there was a “change in the flow characteristics” at the aforementioned threshold wind speed, possibly because continuous wave breaking causes the airflow to separate, effectively limiting the aerodynamic roughness. Although the value of the threshold wind speed at which this regime change occurs is still uncertain (French et al. 2007), there is no evidence to support a continued monotonic increase in CD with wind speed (Bell et al. 2012).
2) Moist enthalpy (sensible and latent heat) coefficients
The exchange coefficient of moist enthalpy Ck (or the coefficients of sensible and latent heat CH and CQ) over the ocean at high wind speeds continues to be problematic. Some of the problem has to do with widespread misconceptions, as Andreas (2011) argues quite convincingly. According to Andreas, there are two pathways for moist enthalpy to be transferred between the air and sea at high wind speeds: the “interfacial” route, involving molecular processes, and the “spray” route, involving transfer from the cooling and evaporation of airborne spray droplets. He also demonstrates that when interfacial transfer alone is considered, CH and CQ are essentially equal to each other and are somewhat weak functions of wind speed; when spray effects are included, the exchange coefficients become strong functions of both wind speed and the air–sea potential temperature difference. Published measurements of CH and CQ appear to agree with the Andreas (2011) calculations for interfacial transfer only, but exhibit large spread at high wind speeds (Drennan et al. 2007; Zhang et al. 2008), which is possibly a manifestation of the spray effect (Bell et al. 2012) or a consequence of sampling error (there is inherent noise in eddy-covariance flux measurements because they try to determine a higher-order moment).
c. Background on flux options in WRF
Because of these uncertainties in CD, CH, and CQ, WRF-ARW allows users to change—through the namelist option isftcflx—how the air–sea fluxes are formulated. The exchange coefficients for drag, sensible heat, and latent heat under neutral stability conditions are, respectively [see (A12)–(A14)],
where k is the von Kármán constant, zref is the reference height (often 10 m), and the subscript N is included as a reminder that these equations are valid only for neutrally stable conditions. It is crucial to note that the surface flux options in WRF-ARW (isftcflx) only affect (5)–(7) by changing the roughness (or scaling) lengths of momentum, sensible heat, and latent heat (z0, zT, and zQ, respectively). This research considers a number of flux options, including all of those that are publicly available in version 3.4.1 of WRF-ARW, as well as some others that are not.
1) Isftcflx options available in version 3.4.1 of WRF-ARW
There are three isftcflx options available in version 3.4.1 of WRF-ARW. For simplicity, each option is referred to by its corresponding number (0, 1, or 2) in the WRF namelist file.
(i) Option 0
Option 0 (hereafter Opt 0) is the default option (and was, until version 3.0, the only option for calculating ocean surface fluxes) in WRF-ARW. The momentum roughness length (all roughness lengths are in meters) is based upon Charnock (1955), and [as a consequence of (5)] yields a monotonic increase in CD with wind speed:
where is the frictional velocity and g the acceleration due to gravity. In the case of no wind ( = 0), the second term in (8) prevents z0 from going to zero. Note that (8) has remained unchanged since 2006 [i.e., several versions of WRF-ARW; Davis et al. (2008)].
For this option, the thermal roughness length zT is set equal to z0 but with the additional constraint that ln(zref/z0) − ψh(zref/L0) ≥ 2 (ψh and L0 are discussed in the appendix). In the case of neutral stability (ψh = 0) and a reference height of 10 m, the constraint simplifies to z0 ≤ 1.35 m, which is satisfied for 10-m wind speeds up to 134 m s−1 (Fig. 1a). It follows from (5) and (6) that CD,N = CH,N for this option (Figs. 1b,c). Additionally, only Opt 0 ignores the effect of dissipative heating on sensible heat flux, which can be significant at high wind speeds (see section 5).
The moisture roughness length zQ takes on a somewhat different functional form:
with the second term including a background molecular diffusivity (Ka = 2.4 × 10−5 m2 s−1). Because the second term in (9) is never negative, zQ ≤ zT = z0 (Fig. 1a) and, thus, from (7) CQ,N ≤ CH,N = CD,N (Fig. 1b).
At this point, it is important to address the peculiarities of Opt 0. Being the first (and default) option, it originated from the MM5 [see section 5.4.3 of Grell et al. (1995)]. Because MM5 was designed for more than just TC prediction, it is unsurprising that its only consideration of oceanic surface fluxes was to calculate z0 by a formula similar to (8) instead of from a land-use category. Consequently, MM5 uses (9) to calculate zQ over both land and water, despite the fact that (9) was developed for an urban–rural canopy (Carlson and Boland 1978). Finally, MM5 (and thus Opt 0) calculates zT (8) differently than zQ (9) because there is molecular diffusion of water vapor, which acts to decrease surface moisture fluxes (Oncley and Dudhia 1995). In our view, the roughness length formulas in Opt 0 are demonstrably inconsistent with a substantial body of research; evidently, this view was shared by a sufficiently large (and/or vocal) group to warrant the development of Options 1 and 2 (hereafter Opt 1 and Opt 2).
(ii) Option 1
This option was developed in response to the finding that CD seemed to level off at hurricane force wind speeds (e.g., Powell et al. 2003; Donelan et al. 2004), and was first implemented in version 3.0 of WRF-ARW. In version 3.4.1, the momentum roughness length for this option is given by a “blend” of the roughness lengths used in versions 3.3.1 and 3.4.0 (R. Torn 2012, personal communication):
Using zQ = 0.95 × 10−4 m from Large and Pond (1982) [note that Davis et al. (2008) erroneously cite Large and Pond (1981)] as a baseline, Opt 1 sets zT = zQ = 10−4 m for all wind speeds. Nevertheless, because CH and CQ are functions of CD [see (6) and (7)], the heat exchange coefficients do vary with wind speed (Fig. 1b).
(iii) Option 2
This is the newest option available and was introduced in version 3.2. The formulation for z0 is given by (10), meaning that Opts 1 and 2 have identical momentum roughness lengths (and thus drag coefficients). The formulas for zT and zQ can be traced back to Brutsaert (1975), and may be expressed as
2) Experimental isftcflx options
Because of the strong sensitivity of simulated TCs to the drag coefficient (e.g., Bao et al. 2012), we ran two additional experiments with different formulas for z0. In both of these experiments, we used (11) and (12) to calculate zT and zQ (i.e., the formulas of Opt 2). We stress that changes to z0 impact the values of zT and zQ (Fig. 1a) even though their formulas are unchanged between these experiments; even so, the differences in CH (and, similarly, in CQ) between these new experiments and Opt 2 are minor, even at high wind speeds.
The first new option—hereafter called PSU—is an ad hoc attempt2 to use the momentum roughness length formula for Opts 1 and 2 at low wind speeds but to allow for a continued increase in drag at wind speeds above ~33 m s−1, yielding
where z0 is given by (10), except without the upper bound of 2.85 × 10−3 m. The value of implementing this option is that the effects of drag at high wind speeds can be isolated (when compared with Opt 2).
The second new option (not shown in Fig. 1) is called K2D0 because it simply uses the thermal roughness lengths from Opt 2 and the momentum roughness length (drag) from Opt 0 [as given by (8)]. K2D0 was tested just for the Katrina case study; results from this option are only discussed within the context of sensitivity to CD. Together, all of these options (summarized in Table 1) provide the most comprehensive investigation of the effects of wind speed–dependent exchange coefficients on simulated TCs to be reported in the literature.
Although the flux options tested are based on WRF-ARW version 3.4.1, we actually ran the simulations using version 188.8.131.52 There were three domains—D01, D02, and D03—with horizontal grid spacings of 27, 9, and 3 km, respectively (all domains used 43 vertical levels); the corresponding time steps were 60, 20, and 20/3 s. Only results from D03 are presented. The present work is mainly focused on Hurricane Katrina, an infamous major hurricane that traversed the Gulf of Mexico in late August 2005 before devastating the northern Gulf coast (Knabb et al. 2006). There are two main components to the experimental procedure: an assimilation stage and a sensitivity stage. These stages are discussed below.
a. Spinup of numerically simulated tropical cyclones
The Katrina simulations were initialized at 0000 UTC 25 August 2005 from Global Forecast System (GFS) initial conditions (ICs) and lateral boundary conditions (LBCs). At this time, an ensemble of 60 forecast members was created by adding random perturbations to the ICs and LBCs; these 60 ensemble members were integrated forward until 1430 UTC to generate a flow-dependent covariance matrix. Between 1430 and 2000 UTC, six rounds of airborne Doppler radar velocity data were assimilated using the ensemble Kalman filter (EnKF) data assimilation technique (Zhang et al. 2009; Weng and Zhang 2012). The ensemble mean at 2000 UTC was integrated forward an additional 4 h to 0000 UTC 26 August 2005. The use of data assimilation required that all grids be fixed in space throughout this period.
Through this point, the experimental procedure is nearly identical to that of Weng and Zhang (2012), with the following exceptions: the use of a newer version of WRF-ARW, a larger horizontal grid in D03 (562 × 562 points), higher resolution (both horizontal and vertical), and isftcflx Opt 2 (as opposed to Opt 0).4 It is crucial to note that the sensitivity to flux option was not tested during the assimilation stage: doing so would have introduced additional uncertainty associated with the treatment of model error by the EnKF algorithm, and accounting for this is beyond the scope of the present paper.
The Grell–Devenyi cumulus scheme (Grell and Devenyi 2002) was implemented for D01, whereas D02 and D03 only used explicitly resolved convection. Also used were the WRF single-moment 6-class (with graupel) microphysics scheme (Hong and Lim 2006), the Rapid Radiative Transfer Model for longwave radiation (Mlawer et al. 1997), and the Dudhia (1989) shortwave radiation scheme. The Yonsei University (YSU) planetary boundary layer scheme (Hong et al. 2006) was used with “MM5 similarity” (sf_sfclay_physics option 1) and five-layer thermal diffusion over land. Single-column ocean mixing was turned off (to isolate the effects of surface fluxes), although we hope to incorporate this important physical process (which would almost certainly result in less intense TCs) at some point in the future.
b. Sensitivity to various flux options
The sensitivity of simulated TCs to various flux options is tested in the second stage of the experiment. These simulations start at 0000 UTC 26 August from the 4-h deterministic forecast based on the EnKF ensemble mean (discussed above). From this point on, the inner two domains were allowed to automatically follow the TC vortex. Five physics sensitivity experiments5—one for each of the flux configurations described in section 2—were integrated forward 120 h to 0000 UTC 31 August. All results presented herein are from within this 5-day window.
4. Temporal evolution of simulations (12–72 h)
The tracks and common metrics of TC intensity—minimum SLP and maximum 10-m wind speed—over the entire 120-h sensitivity stage are plotted in Figs. 2a–c using the output of the WRF automatic vortex-following algorithm. All simulations moved slower than, and to the left of, the observed best track. Generally, the simulations continued to intensify up until landfall whereas Katrina was at its peak ~18 h before landfall; this discrepancy may be partially due to the neglect of ocean feedbacks in the model. Regardless, the assimilation of airborne Doppler radar data yields remarkably accurate forecasts of track, SLP, and wind speed (Weng and Zhang 2012). With the exception of these three metrics, data analysis was limited to hours 12–72 of the sensitivity experiments (1200 UTC 26 August–0000 UTC 29 August) to exclude times when the simulated storms were near land.
a. Intensity (wind speed and SLP)
Figures 2b and 2c demonstrate the strong sensitivity of intensity to the surface flux option. Opt 1 is by far the most intense in terms of both wind speed and SLP, which is to be expected because, at high wind speeds, it has the largest Ck/CD ratio (Fig. 1c) and shares (with Opt 2) the lowest CD (Fig. 1b). Unsurprisingly, the options in which CD does not level off have the weakest wind speeds. The minimum SLP is negatively correlated with CQ/CD; however, maximum wind speed is not positively correlated with CQ/CD. Thus, the pressure–wind relationship varies between flux options (as will be discussed in detail later on).
b. Track, size, and strength (12–72 h)
Track is not very sensitive to the flux option (Fig. 2a), but it appears that stronger storms (Opt 1) fall on the left side of the track envelope whereas weaker storms (PSU) fall on the right side. Because the tracks are all very similar, the differences in underlying sea surface temperature (which remains fixed throughout the simulation because there was no ocean feedback) between the options are negligible. Consequently, differences in TC intensity and structure are solely due to the differences in surface flux formulas.
Size and strength—defined as the radius of tropical storm force winds [34 kt (~17.5 m s−1)] and the relative angular momentum averaged inside a radius of 300 km, respectively (Holland and Merrill 1984)—both increase linearly with time between hours 12 and 60 (Figs. 2d,e). But after 60 h (1200 UTC 28 August), size appears to grow more slowly whereas strength continues to increase. The intensity metrics (except for Opt 1) level off by ~60–66 h (Figs. 2b,c), which might suggest that the simulated TCs reach a quasi-steady state. Nevertheless, our use of a real data case (as opposed to an idealized simulation) makes us hesitant to claim that the TCs are in a quasi-steady state for the purpose of supporting or refuting PI theory [see Bryan (2013) for an excellent discussion of this topic]. Additionally, despite significant intensity differences among the various simulations, the size and strength parameters are more or less the same (except for Opt 1), which hints at fundamental differences in the radial distributions of near-surface winds.
c. Radial variability (12–72 h)
Radius–time Hovmöller diagrams of azimuthally averaged fields allow for the investigation of radial variability (i.e., structural differences) both within and between simulated TCs during the intensification period (between 1200 UTC 26 August and 0000 UTC 29 August 2005). Of particular interest are the low-level wind fields, moist enthalpy fluxes, and SLP.
1) Low-level winds
Because Figs. 2d and 2e suggest that there are fundamental differences in the distributions of low-level winds between flux options, radius–time plots of 10-m tangential and radial winds are given in Fig. 3; these plots confirm that the wind field expands throughout the analysis period. Two interesting points emerge from Fig. 3. First, the inflow of Opt 0 is closer to the TC center than in the other simulations, especially by 0000 UTC 29 August (72 h). Second, the inflows of Opt 2 and PSU are somewhat similar. This is because drag directly causes inflow, whereas moist enthalpy fluxes are a secondary effect (higher moist enthalpy fluxes yield more intense storms, which yield stronger inflow).
Differences in the 10-m wind speed between Opt 2 and each of the three other extensively tested options (0, 1, and PSU) are shown in Figs. 4a–c. Between Opts 1 and 2 (Fig. 4a), significant differences (greater than 2 m s−1) extend out to almost 150 km by 72 h, which is not surprising because these two options have the same CD but different CH and CQ. The wind speed differences between PSU and Opt 2 (Fig. 4c) are confined mainly inside ~90 km, which is to be expected because these options have nearly identical values of CH and CQ (Fig. 1b).
The 10-m zonal wind is diagnosed in WRF from the prognostic winds at the lowest model level using the following formula:
where subscripts 10 and LL denote reference heights of 10 m and the lowest model level, respectively; u is the zonal wind; and CD is evaluated from (A9) at the appropriate reference height. Note that (14) can also be applied to the meridional wind. To ensure that our results are not just a consequence of the diagnosis, we made plots similar to Figs. 3 and 4a–c but for winds at the lowest model level (not shown) and found that the various flux options affect the prognostic winds in a manner similar to the 10-m winds.
2) Moist enthalpy fluxes and SLP
Radius–time plots (Fig. 5) show that sensible heat fluxes are smaller than latent heat fluxes, a result that is consistent with those derived from field measurements (e.g., Zhang et al. 2008). The relationship between SLP and surface heat fluxes (greater fluxes mean more warming, which yields lower pressure) is quite evident in our results (Figs. 4d–f and Fig. 5). In fact, the effect of heat fluxes on pressure extends well beyond the TC center, as will become more apparent in the next section.
5. Results at 1200 UTC 27 August and 0000 UTC 29 August 2005
Even though the simulated TCs might not reach the quasi-steady state necessary for PI considerations (Bryan 2013), focusing on single times can still be beneficial. For this purpose, 1200 UTC 27 August and 0000 UTC 29 August—36 and 72 h into the sensitivity stage, respectively—were chosen. At 36 h, the simulated TCs were at their southernmost location and undergoing steady intensification, whereas at 72 h the TCs were south of the Louisiana coast and only Opt 1 was still intensifying (Figs. 2a–c). As in section 4, all results are azimuthal averages.
a. Radial distributions of 10-m wind speed, surface fluxes, and SLP
The radial distributions of 10-m wind speed, surface heat fluxes, SLP, and its radial gradient for these two times are plotted in Fig. 6. Note that this figure also includes K2D0, which has the same z0 as Opt 0 and the same formulas for zT and zQ as Opt 2 and PSU. The radius of maximum wind at 36 h (Fig. 6a) ranges between 24 and 30 km for the five experiments, but by 72 h a clearer picture emerges (Fig. 6b) as to how drag and enthalpy affect low-level wind speed. Opt 0 and K2D0 are peaked at a radius of 21 km, whereas Opts 1 and 2 are peaked at 33 and 30 km, respectively. The PSU option, which at high wind speeds has values of CD between Opts 0/K2D0 and Opts 1/2, is peaked at 24 km. These results seem to suggest that CD can modulate the radius of the strongest 10-m winds.
The intensity of the low-level wind speed obviously governs the sensible and latent heat fluxes (Figs. 6c,d). The fluxes are generally related to CH and CQ, except that CH does not account for dissipative heating, which is included in all options6 except for Opt 0. Figures 6c and 6d also plot what the sensible heat flux for Opt 0 would have been if the dissipative heating term were included.7 With this change, the relationship between exchange coefficients and heat fluxes is perfectly straightforward: higher coefficients yield higher fluxes for a given wind speed [e.g., (A2) and (A3)].
As previously noted, heat fluxes affect the pressure distribution beyond the inner core of the TC. At large radii, the options with the lowest pressure, Opts 0 and 1 (Figs. 6e,f), also have the highest heat fluxes (Figs. 6c,d). By 72 h (Fig. 6f), the situation gets more complicated inside 50 km because Opt 0 and K2D0 have narrow but sharp peaks in the radial pressure gradient; thus, in these simulations of Katrina, CD can affect the shape of the pressure gradient field and the size of the inner core (as measured by the radius of maximum 10-m winds).
The combined effects of CD on low-level winds and minimum SLP have major implications for the relationship between these two intensity metrics. For example, at 72 h, K2D0 has lower SLP than PSU only inside of 15 km from the TC center, but PSU has a slightly higher maximum 10-m wind speed (53.2 versus 50.9 m s−1). Thus, surface drag affects the pressure–wind relationship of a relatively mature TC and shifts the radius of maximum 10-m winds.
b. Radius–height structure
Radius–height plots8 of the tangential and radial winds at 72 h (Fig. 7) indicate that the different flux options yield the same general structure, albeit with some differences. These differences are better seen in Fig. 8, which zooms in to the lowest 4 km above the surface and only 150 km outward from the TC center. Looking at Fig. 8, PSU has the weakest winds (tangential winds barely reaching 70 m s−1) whereas Opt 1 has the strongest winds (tangential winds greater than 90 m s−1). Additionally, the options with higher CD, Opt 0 and PSU (Figs. 8c,d), have a secondary maximum of inflow close to the eyewall (~25 km). The stated findings also hold at 60 h (1200 UTC 28 August; not shown), which means that the differences in Figs. 7 and 8 cannot be fully explained by the TCs being at various stages of an eyewall replacement cycle.
Clearly, the choice of surface flux option has a tremendous impact on simulations of Hurricane Katrina. It is reasonable, then, to ask how applicable the results from this one case are to TCs in general. To answer this question, we tested Opts 0, 1, 2, and PSU (as formulated in version 3.3.1 of WRF-ARW)9 on deterministic forecasts of nearly all North Atlantic TCs that were sampled by airborne Doppler radars during the 2008–11 hurricane seasons (for a total of 71 cases). As with the Katrina simulations, these forecasts were run after a period of assimilating airborne Doppler radar data by the EnKF algorithm (flux option was not tested during the data assimilation stage). The resulting pressure–wind relationships from these experiments, when aggregated to combine all forecasted lead times for all simulations, are shown in Fig. 9 along with the pressure–wind relationships for the observed best tracks of these TCs. Despite considerable spread in the data, the lines of best fit bolster our findings from Katrina. First, in agreement with recent work (Bao et al. 2012), CD appears to affect the pressure–wind relationship: for a given pressure, options with higher drag (such as Opt 0) will have weaker 10-m wind speeds. The moist enthalpy fluxes, on the other hand, have less of an impact on the pressure–wind relationship, but still control both metrics of TC intensity: greater fluxes yield storms with faster winds and lower central pressures (cf. Opts 1 and 2).
As stated above, Fig. 9 provides strong evidence to support the notion that changing the surface flux option systematically impacts the intensity of simulated TCs. Nevertheless, it is necessary to ensure that the differences in intensity (at 72 h, say) are not solely a consequence of the projection of slight differences in state variables (which arise in the first few hours of forward integration) onto growing error modes. To do this, we randomly chose 10 of the 60 ensemble forecast members at 2000 UTC 25 August (immediately following the final EnKF update) and reduced the differences from the ensemble mean by 90% in D03 (the ICs in the outer domains were set to those of the ensemble mean). These 10 members were integrated forward using Option 2 until 0000 UTC 31 August. The results (not shown) are unsurprising: Opts 0 and 1 have significantly different intensities than any of the perturbed members, whereas PSU has significant differences only for very strong wind speeds (i.e., the regime in which CD for PSU is not equal to CD of Opt 2; see Fig. 1). In addition, the pressure–wind relationships of the 10 IC sensitivity experiments, while very close to one another and to Opt 2, are substantially different than the pressure–wind relationships of Opts 0, 1, and PSU.
Another important question is which flux option is the “best.” Before going any further, we cannot stress enough that our simulations were unable to incorporate any ocean coupling (one- or three-dimensional); accounting for TC-induced ocean cooling will almost certainly reduce intensity forecasts across all options. That said, the answer obviously depends on what the WRF user is trying to accomplish. Forecasters in an operational setting are most concerned with 10-m wind speeds because this metric defines the Saffir–Simpson hurricane wind scale, which is relied upon heavily by government agencies and the general public. In terms of root-mean-squared errors and biases, the best options for 10-m wind speed forecasts are Opt 0 and PSU (Y. Weng 2012, personal communication).10 On the other hand, those in a research setting may be interested in a better pressure–wind relationship because it would suggest a more physically realistic TC; in this case, Opts 1 and 2 are favored. Regardless, these recommendations may fail if other components of the numerical model (such as planetary boundary layer parameterization, horizontal grid spacing, and particularly the inclusion of ocean coupling) are different from the experimental setup discussed in section 3. We tested the Mellor–Yamada Nakanishi Niino (MYNN; Nakanishi and Niino 2006) planetary boundary layer scheme with Opts 0, 1, and 2 for Katrina (not shown), and found the pressure–wind relationships were very similar to the YSU simulations. Obviously, more rigorous testing is required before our results can be generalized across different boundary layer schemes.
The accuracy of TC intensity forecasts has stagnated over the past two decades, despite significant progress in the skill of TC track forecasts. Previous research involving a combination of theory, observations, and numerical modeling has found that TC intensity is quite sensitive to fluxes of momentum and moist enthalpy across the air–sea interface. Unfortunately, NWP models cannot resolve these fluxes directly and thus must resort to imperfect parameterization schemes.
This study used WRF-ARW to examine the impacts of different surface flux parameterizations (including some that are not publicly available) on the intensity and structure of Hurricane Katrina (2005). Although TC track was not significantly affected, the common intensity metrics—minimum SLP and maximum 10-m wind speed—varied widely among the flux options, which is consistent with the results of others. As expected (e.g., Ooyama 1969; Rosenthal 1971; Bao et al. 2002), options with higher moist enthalpy fluxes yielded storms with deeper pressures and stronger winds. The contribution of dissipative heating to surface sensible heat flux is not negligible when the near-surface winds are at hurricane strength.
In agreement with recent work by Bao et al. (2012), the pressure–wind relationship was affected by the drag coefficient CD. Specifically, we found that capping CD at a constant value for very high wind speeds, as suggested by observational data (Powell et al. 2003; Donelan et al. 2004), yields a pressure–wind relationship quite similar to that of the observed best track; this result still holds when the dataset was expanded to include several North Atlantic TCs from the 2008 to 2011 seasons. The more CD was allowed to increase at high wind speeds, the worse the pressure–wind relationship became. The most likely physical explanation is that for given pressure and low-level wind fields, higher drag weakens the 10-m wind field as a direct consequence of (14).
The TC structure was also impacted by the choice of the surface flux scheme. For example, an increase in the exchange coefficients of sensible heat CH and latent heat CQ (and thus the fluxes themselves) resulted in a larger storm with higher relative angular momentum. Changes to CD, on the other hand, mainly impacted the TC within ~90 km of the center. For these simulations of Katrina, increased drag yielded a sharper radial pressure gradient, a smaller radius of 10-m maximum winds, and a secondary low-level inflow maximum.
Our research does not solve the ongoing problem of accurately forecasting TC intensity; rather, it serves to raise awareness of how (and, to some extent, why) different ways of parameterizing fluxes across the air–sea interface can affect simulated TCs. Furthermore, there are other parameterized processes that have a significant influence on TCs (e.g., Bao et al. 2012; Bryan and Rotunno 2009; Bryan 2012; Kepert 2012). Continued improvements in high performance computing will allow numerical models to have even smaller grid spacing, to the point where turbulent large eddies will be explicitly resolved (Rotunno et al. 2009). Even with large eddy simulations, accurate representation of fluxes will remain a challenge because of the complexity of surface waves and associated spray—problems that will require coupling with wave and three-dimensional ocean models. Such continuous change means that attempting to find a set of physical parameterization schemes that is optimal for TC prediction would be a fruitless endeavor.
One possible way to improve a physics scheme without worrying about changes to other parts of the model configuration is the implementation of parameter estimation, which is a data assimilation technique that estimates model parameters in addition to atmospheric state variables (Anderson 2001). Parameter estimation has shown promising results in tests of both simple models (Aksoy et al. 2005, 2006a) and full-physics models (Aksoy et al. 2006b; Hu et al. 2010). The next step of this research will incorporate parameter estimation into the WRF surface layer code, with a focus on the flux options, to determine if TC intensity forecasts can be improved.
This material is based upon work supported by a National Science Foundation (NSF) Graduate Research Fellowship under Grant DGE-0750756. Additional support is provided by ONR Grant N000140910526, NSF Grant ATM-0840651, and NOAA (HFIP). We appreciate the valuable comments provided by George Bryan, Kerry Emanuel, and two anonymous reviewers on earlier versions of the manuscript. The authors acknowledge the Texas Advanced Computing Center (TACC) at the University of Texas at Austin (http://www.tacc.utexas.edu) for providing high performance computing resources that have contributed to the research results reported within this paper. Yonghui Weng is thanked for providing the data used to create Fig. 9, and Ryan Torn is thanked for sharing the surface layer code of version 3.4.1 before its public release. Discussions with Daniel Stern were extremely beneficial, as were those with Chris Davis, Jimy Dudhia, and Wei Wang.
Atmospheric Surface Layer
The fluxes of momentum τ, sensible heat H, and latent heat E across the atmospheric surface layer are
where ρ is the density of air; is the friction velocity; and are the surface layer temperature and moisture scales, respectively; Δ(U, θ, q) are the respective differences in wind speed, temperature, and water vapor between a reference height zref (often 10 m) and the bottom of the surface layer (note that U = 0 at the bottom of the surface layer); cp is the specific heat capacity of air; Lυ is the enthalpy of vaporization; and CD, CH, and CQ are the respective bulk exchange coefficients for drag, sensible heat, and latent heat.
Monin–Obukhov similarity theory may be used to calculate , , and :
where k is the von Kármán constant; z0, zT, and zQ are the roughness lengths for momentum, sensible heat, and water vapor (latent heat), respectively [see (8)–(13) in main text]; ψm and ψh are the stability correction functions for momentum and heat; and L0 is the Obukhov length, given by
where g is the acceleration due to gravity and θ0 is the base-state temperature. In neutral stability,
b. Exchange coefficients
where the subscript N is included as a reminder that these equations are valid only for neutrally stable conditions.
Zhang et al. (2011) obtained very good forecasts of 10-m wind speed (but not of pressure) using Opt 0; offline tests of Opts 1 and 2 improved the pressure–wind relationship at the cost of less accurate wind speed forecasts. The PSU option was developed to try and draw from the strengths of Opt 0 and Opts 1/2.
We ran the simulations using V3.4.0 because V3.4.1 was made publicly available very late into this research. Fortunately, Ryan Torn (2012, personal communication) provided us beforehand with the surface layer code that was implemented in V3.4.1.
Actually, the assimilation stage used Opt 2 as it was formulated in V3.4.0; in this earlier version, z0 is given by (10c) and retains the upper bound of 2.85 × 10−3 m. Time and computational constraints prevented us from rerunning the entire assimilation stage with Opt 2 as formulated in V3.4.1, although we believe that this difference is not important because physics sensitivity was tested after the assimilation period.
Dissipative heating was not included in Opt 1 before the release of version 3.2.
Obviously, had dissipative heating been included in Opt 0 from the beginning (i.e., during model integration), the nonlinear interaction between sensible heat fluxes and the atmospheric state would have produced a significantly different TC and thus different heat fluxes by 36 h.
For these plots, data were interpolated to a vertical grid spacing of 100 m using a piecewise cubic Hermite polynomial.
This older version was used so as to be consistent with ongoing research (Y. Weng 2012, personal communication) that adds 12 new cases from 2011 (7 from Hurricane Irene and 5 from Hurricane Rina) to the 59 cases examined in Zhang et al. (2011). In V3.3.1, Opts 1, 2, and PSU use (10d) to calculate z0. The lower bound of 1.27 × 10−7 is used for all three options; Opts 1 and 2 use an upper bound of 2.85 × 10−3 whereas PSU still uses (13) to remove the upper bound. We will test the V3.4.1 formulations in the near future.
These simulations did not include single-column ocean mixing.