## Abstract

A synthetic dataset is used to show that apparent variations between different stability classes in the mean drag coefficient, *C*_{D10n}, to wind speed relationship can be explained by random errors in determining the friction velocity *u*∗. Where the latter has been obtained by the inertial dissipation method, the variations in *C*_{D10n} have previously been ascribed to an imbalance between production and dissipation in the turbulent kinetic energy budget. It follows that the application of “imbalance corrections” when calculating *u*∗ is incorrect and will cause a positive bias in *C*_{D10n}, by about 10^{−4}.

With no imbalance correction, random errors in *u*∗ result in scatter in the *C*_{D10n} values, but for most wind speeds, there is no mean bias. However, in light winds under unstable conditions random errors in *u*∗ act to positively bias the calculated *C*_{D10n} values. This is because the stability related effects are nonlinear and also because for some records for which *C*_{D10n} would be decreased, the iteration scheme does not converge. The threshold wind speed is typically 7 m s^{−1}, less for cleaner datasets. The biased *C*_{D10n} values can be avoided by using a *u*∗ value calculated from a mean *C*_{D10n}–*U*_{10n} relationship to determine the stability. The choice of the particular relationship is not critical. Recalculating previously published *C*_{D10n} values without imbalance correction, but with anemometer response correction, results in a decrease of *C*_{D10n} but only by about 0.05 × 10^{−3}.

In addition to removing a previous cause of scatter and uncertainty in inertial dissipation data, the results suggest that spurious stability effects and low wind speed biases may be present in *C*_{D10n} estimates obtained by other methods.

## 1. Introduction—The apparent need for an imbalance term

Over a number of years we have been using the inertial dissipation method to determine the wind stress over the ocean [e.g., Yelland and Taylor 1996; Yelland et al. 1998, henceforth YT96 and Y98, respectively]. The aim has been to determine the 10-m neutral values for the drag coefficient *C*_{D10n}, as a function of the wind speed *U*_{10n}, under a wide range of conditions:

with the friction velocity *u*∗ determined by the inertial dissipation technique (e.g., Fairall and Larsen 1986):

where *U*_{rel} is the relative wind speed measured at the anemometer site (at height *z* above the sea), and *k* and *K* are the von Kármán and Kolmogorov constants, respectively. The power spectral density (PSD) of the wind component *S*_{u}(*f*) is determined at a spectral frequency, *f,* which is typically >2 Hz. Since this frequency is above those induced by ship motion, the method is particularly attractive for use at sea. However, the method has been criticized (e.g., Donelan et al. 1997) because of the assumptions needed with regard to the form of the nondimensional function *ϕ*_{ɛ} for the dissipation, ɛ:

where *ϕ*_{ɛ} is determined from the nondimensionalized turbulent kinetic energy (TKE) balance:

The terms on the right-hand side of Eq. (4) respectively represent the mechanical production, the buoyant production, and the vertical divergence of the turbulent transport and pressure transport terms; all are expected to be functions of the stability parameter *ζ*:

where the Obukhov length, *L,* is a function both of *u*∗ and of the vertical flux (*F*_{Tυ}) of virtual temperature *T*_{υ}:

There are various difficulties in solving Eqs. (1)–(6) to obtain *u*∗. First, a direct determination of *F*_{Tυ} often is not available and the flux is calculated from the bulk formulas for the fluxes of sensible and latent heat. Second, since *ζ* depends on *u*∗ (and the bulk formulas depend on *ζ*) iteration is required. This is typically performed by assuming neutral stability and reiterating for *u*∗, *U*_{10n}, and *ζ*; however, the iteration may not always converge (section 2c). Third, the exact forms of *ϕ*_{m}(*ζ*), *ϕ*_{t}(*ζ*), and *ϕ*_{p}(*ζ*) are open to debate. Following the review of Edson et al. (1991), YT96 and Y98 used

However, more recently, Edson and Fairall (1998) have suggested different forms for *ϕ*_{m}(*ζ*), which have the correct form in the free convection limit.

With regard to the vertical divergence terms, Large (1979) argued that the available evidence suggested that

which implies a balance between dissipation and the sum of mechanical and buoyant production:

However, using this assumption in our analysis resulted in values for the mean *C*_{D10n} to *U*_{10n} relationship, which apparently varied with stability. This led to the introduction (YT96) of an “imbalance term,” *ϕ*_{D}, which was apparently a function both of *ζ* and wind speed. The form of the imbalance term was determined empirically so that the mean *C*_{D10n}–*U*_{10n} relationships for different ranges of *ζ* were all collapsed as far as possible on to the relationship obtained from near-neutral data. Using a different implementation of the inertial dissipation method, and data from a different experiment, Dupuis et al. (1997, henceforth DTWK97) similarly identified a need for an imbalance term. Following previous studies of the TKE budget (e.g., Thiermann and Grassl 1992; Vogel and Frenzen 1992; Wyngaard and Cote, 1971), both YT96 and DTWK97 assumed that this term represented an imbalance between production and dissipation of TKE. Thus

The magnitude of the required imbalance term was found to be significantly larger than the uncertainty in the formula for *ϕ*_{m}(*ζ*). Recently Edson and Fairall (1998) have suggested new forms for the various terms in the TKE budget also including a significant imbalance term.

Unfortunately, analysis of further datasets has shown that the imbalance term proposed (e.g., in YT96) does not remove the stability dependence of the results. In the remainder of this note we will construct a synthetic dataset and show that the imbalance term arises from the dependence of *ζ* on *u*∗. Furthermore, that random variations in *u*∗, as might be caused by measurement errors, result in a form of *ϕ*_{D} similar to that suggested in previous studies. We will then suggest an implementation of the inertial dissipation method, which avoids the problem of the imbalance term and correctly processes the synthetic data. Finally we will test the new program using the YT96 data and determine the difference from the *C*_{D10n}–*U*_{10n} relationship published by Y98 using the same dataset.

## 2. Examination of the problem using the synthetic data

### a. Construction of the synthetic dataset

We shall make use of data obtained during the HEXMAX experiment (Smith et al. 1992) on the R/V *Frederick Russell.* Mean values of air and sea temperature, humidity, wind velocity, and pressure were merged with inertial dissipation data, which had been obtained by researchers from the Naval Postgraduate School (NPS) using sensors mounted on a 10-m mast erected on the ship’s foredeck. For constructing the synthetic dataset we will ignore the NPS data and use only the mean meteorological values as a representative dataset collected at sea. For our present purposes any other marine dataset could be used to obtain similar results. However, the NPS data will later be of value (section 3b).

The synthetic data were constructed by using a standard bulk formula program (as used by Josey et al. 1999) to calculate the *U*_{10n} and *u*∗ values corresponding to the mean values of true wind, temperatures, humidity, and pressure data. The Smith (1980) drag coefficient and Smith (1988) transfer coefficients for sensible and latent heat, *C*_{T} and *C*_{Q}, were used:

From the calculated *U*_{10n} and *u*∗, values for *ϕ*_{m}(*ζ*) and *S*_{u}(*f*) were calculated [Eqs. (2) and (6)]. Thus these *S*_{u}(*f*) values represent those that would be measured in the absence of any measurement error if the Smith (1980) formulas represent the exact values of the drag coefficient and the Smith (1988) heat and water vapor transfer coefficients are the correct values for the heat fluxes. These will be referred to as the synthetic “perfect” data.

A set of simulated “observed” *S*_{u}(*f*) values was also created by adding white noise “error” values [randomly distributed between ±20% × *S*_{u}(*f*)] to the perfect *S*_{u}(*f*) values, the resulting rms scatter was 15%. For comparison, the results of Yelland et al. (1994), who compared PSD values from four anemometers on the same mast, suggest an rms scatter of about 10% for sonic anemometer data and 20% for a propeller/vane-type instrument. Thus for sonic anemometer data the assumed random error may overestimate the magnitude of the instrumental measurement error in *S*_{u}(*f*). However, there is evidence (see section 4 below) that variations in airflow distortion by the ship due to changing relative wind direction (Y98) will increase the random errors and we have also assumed that the height of measurement, true wind speed, air temperature, humidity, sea temperature, etc., are all exactly known. The effect of various measurement errors on the determination of *C*_{D10n} using the inertial dissipation method was discussed in more detail by Yelland et al. (1994); here we will simply suggest that the magnitude of the random error assumed is similar to that observed.

### b. Analysis of the synthetic data

First, we applied our inertial dissipation analysis scheme (as described in YT96 but with no imbalance term) to the synthetic perfect data. As might be hoped, all the calculated *C*_{D10n} values exactly corresponded to the Smith (1980) formula [Eq. (11a)]. The iteration converged for all data points and the results did not show any dependence on stability. We then applied the same program to the synthetic “observed” dataset. The *C*_{D10n} to *U*_{10n} plot (Fig. 1) now showed partitioning by stability with, for *ζ* < 0 and higher wind speeds, values of *C*_{D10n} lying below the Smith (1980) relationship. Conversely, at wind speeds below about 6 or 7 m s^{−1}, some *C*_{D10n} values were increased significantly above the Smith (1980) relationship.

### c. Increased *C*_{D10n} values at low wind speeds

Consider first the apparent increase of *C*_{D10n} at low wind speeds. This was found to be due to two effects:the nonlinearity of the *C*_{D10n}–*S*_{u}(*f*) relationship (caused by the apparent changes in stability) and also nonconvergence of the iteration procedure in the program (which occurred for 31 points). The nonlinearity meant that a random increase in *S*_{u}(*f*) increased *C*_{D10n} more than a random decrease. The nonconvergence is illustrated in Fig. 2. The *ζ* values for data points below 5 m s^{−1} are plotted as a function of the random increase or decrease in *S*_{u}(*f*). For all stable values (*ζ* > 0), and for values where the *S*_{u}(*f*) value was increased (leading to higher *C*_{D10n} values), the iteration converged. However, for most of the unstable points for which the *S*_{u}(*f*) value was randomly decreased, the iteration did not converge. The result was an apparent increase of *C*_{D10n} compared to the true value at low wind speeds.

For the unstable data, whereas DTWK97 suggested an imbalance term of −0.65*ζ,* which increased the calculated *C*_{D10n} values, the imbalance term of YT96 [(2 − *U*_{10n}/3)*ζ*] agreed with the DTWK97 value at about 7 m s^{−1} but decreased *C*_{D10n} for wind speeds below 6 m s^{−1}. DTWK97 suggested that the reason was that YT96 had not taken the effects of nonconvergence into account when determining the imbalance term. Indeed, by restricting the accepted relative wind directions to ±10° of the bow, Yelland (1997, henceforth Y97) reduced the noise in the YT96 dataset. The result was fewer cases of nonconvergence and a modified imbalance term, which only reduced the *C*_{D10n} values for wind speeds below about 3 m s^{−1}. However, Y98 restricted their analysis to data at winds speeds over 6 m s^{−1} and used the YT96 imbalance term.

### d. Apparent stability effects at higher wind speeds

Consider now the effect of a random change in *S*_{u}(*f*) in the higher wind speed ranges. The change in *S*_{u}(*f*) causes the calculated *u*∗ value to change. The transfer coefficients (and therefore roughness lengths) for sensible and latent heat are assumed constant. The change in *u*∗ therefore changes the Obukhov length *L* (which varies with *u*^{3}_{∗}). For example, an increase in *u*∗ results in an apparent increase of turbulence over buoyancy forces and a decrease in the magnitude of the calculated *ζ.* But the increased *u*∗ also results in a larger apparent value of *C*_{D10n} (which varies with *u*^{2}_{∗}). Thus randomly increased *C*_{D10n} values are always associated with randomly decreased |*ζ*| and vice versa. Since, at these wind speeds, the relationship between the change in *S*_{u}(*f*) and the calculated *C*_{D10n} is relatively linear, the random increases and decreases of *C*_{D10n} will tend to cancel and the mean value at any particular wind speed should not change. However, because the air–sea temperature difference over the open ocean is typically small, at these higher wind speeds, the distribution of *ζ* is skewed toward neutral data. Thus in any stability range there are always a greater number of data points for which a decrease in *S*_{u}(*f*) (associated with decreased *C*_{D10n}) will result in the calculated *ζ* being moved farther away from neutral (and into a different range) than there are points which will be moved closer to neutral. Thus biases will appear if the data are subsetted by stability as in Fig. 1. This effect is not confined to data obtained by the inertial dissipation method. It will be present in any study in which data depending on *u*∗ are subsetted according to stability.

### e. The simulated imbalance term

The required “imbalance term” implied by our simulated data can be obtained from

where all the terms on the right-hand side are evaluated from the simulated observed data except *u*∗ for which the synthetic perfect value is used. By subsetting the results from the simulated observed data into different wind speed ranges we can compare with the wind speed–dependent imbalance term suggested by YT96. The results are shown in Fig. 3. For unstable data the agreement is good, particularly for the 7–11 m s^{−1} range, which included much of the data. The increasing slope of the simulated imbalance with increasing wind speed is similar to the YT96 imbalance term. For stable data the agreement is poor in the low wind speed range, 3–7 m s^{−1}. However YT96 had little data in this category and were unsure of the required imbalance. Again the agreement is good for the 7–11 m s^{−1} range.

Various published formulas for the imbalance term in near-neutral conditions are shown in Fig. 4 together with the simulated imbalance for wind speeds above 10 m s^{−1} (chosen to remove any possible low wind speed effect). As previously noted, the use of a cleaner dataset resulted in the magnitude of the Y97 imbalance term being less than that suggested by YT96. Most of the simulated data points lie between DTWK97 and the Edson and Fairall (1998) formulas. For *ζ* < −0.1 the simulated imbalance is similar to the published values. However, closer to neutral the simulated imbalance is small. This discrepancy would be greater if simulated values from lower wind speeds were added. Nevertheless, we believe that the comparisons shown in Figs. 3 and 4 are evidence that it is incorrect to use an imbalance term in order to remove the apparent stability dependence. Indeed, we will show that, for wind speeds above about 6 m s^{−1}, it will positively bias the mean *C*_{D10n} values (see section 4).

## 3. A modified form of the inertial dissipation method

### a. Method of stability calculation

While use of the imbalance term may be incorrect, processing simulated “observations” without application of the imbalance term has been shown to result in nonconvergence and artificially increased *C*_{D10n} values at wind speeds below about 6 m s^{−1}. A remedy is to change the way *ζ* is calculated during the iteration. Instead of using the *u*∗ value calculated from *S*_{u}(*f*) to determine *ζ,* the calculated *U*_{10n} value is used in a bulk *C*_{D10n}–*U*_{10n} relationship to determine a “bulk” *u*∗ value. This bulk *u*∗ value is then used in the *ζ* determination and also used to calculate the 10-m neutral values of temperature and humidity. However, while the *ζ* determination is based entirely on bulk formulas, it is the inertial dissipation derived *u*∗ value that is used to calculate the next *U*_{10n} value, as before. A somewhat similar scheme was suggested by Large (1979) and Large and Pond (1981) to avoid noise in the buoyancy flux data.

### b. Results

A new program was created that used the Smith (1980) ,*C*_{D10n}–*U*_{10n} relationship [Eq.(11a)] for calculating the bulk *u*∗ value. Applying our new program to the simulated observed *S*_{u}(*f*) gives the result in Fig. 5; all records converged and, as expected, there is no apparent effect of stability. The enhanced *C*_{D10n} values at low wind speeds have been avoided. Values for one-way regressions calculated from the individual data points are shown in Table 1. For the synthetic perfect data the new program exactly retrieved the Smith (1980) relationship. For the simulated observational data both the original and new programs successfully retrieved the slope and intercept of the Smith relationship within the calculated error and to similar accuracy. However, the scatter was slightly reduced using the new program.

Since the new program used the Smith (1980) formula to calculate stability, did that force the *C*_{D10n} results toward the Smith (1980) relationship? Figure 6 shows the results for *C*_{D10n}–*U*_{10n} calculated using the new program applied to the actual NPS *S*_{u}(*f*) data. Also shown are results from running a different version of the new program that calculated *ζ* using a HEXMAX *C*_{D10n}–*U*_{10n} relationship [Eq. (13) in Smith et al. 1992]. The two sets of results are almost identical and similar to the chosen HEXMAX relationship; this demonstrates that the particular choice of the bulk *C*_{D10n}–*U*_{10n} relationship is shown to be not important.

### c. Convergence

Using the NPS *S*_{u}(*f*) data and our original inertial dissipation program with no imbalance term, 46 out of the total of 2811 records did not converge. Using the Y97 imbalance term reduced the nonconvergence to 15 records. With the new program two or four records did not converge depending on whether the Smith (1980) or Smith et al. (1992) bulk relationship was used. These data were the very lowest wind speed records. Accepting that extrapolation of the Smith (1980) relationship to zero wind speed is not justified, the formulas were modified to give an increase in *C*_{D10n} as *U*_{10n} approached zero:

Using this modification, no records failed to converge.

### d. Implications for the Yelland et al. (1998) *C*_{D10n}–*U*_{10n} relationship

The processing scheme used by Y98 (i.e., with the imbalance term developed by YT96) was applied to the simulated observed data. The results of linear regression on the individual points are shown in Table 1. The calculated *C*_{D10n} values are biased high; the mean bias was 0.08 ± 0.004 (standard error). Thus using the new program, lower values for *C*_{D10n} will be obtained. However, since the work described in Y98 was completed, Henjes et al. (1999) have demonstrated that a correction is needed for the pulse averaging applied within the sonic anemometer. The data used by Y98 was a mixture of recording runs sampled at 56 and 21 Hz. The *C*_{D10n} values calculated from the 56-Hz data need little correction; however, those recorded at 21 Hz are increased by about 10^{−4}. Inspection showed that this difference was apparent in the original datasets, confirming the need for the Henjes et al. (1999) correction. The overall change to the Y98 results by using the new program and the pulse averaging corrections is shown in Fig. 7. The mean *C*_{D10n} values have been decreased by about 0.05 × 10^{−3}. The difference between the new relationship,

and that quoted by Y98 is not significant. Figure 7 also demonstrates that, even with the new program, large *C*_{D10n} values are observed at low wind speeds. Since the new program did not introduce high values using the synthetic data, these high *C*_{D10n} values may be real. However, sources of error not considered here (e.g., the wind speed determination) also have the potential to produce errors in *C*_{D10n}.

## 4. Summary

Because both the drag coefficient *C*_{D10n} and the stability parameter *ζ* are proportional to the friction velocity *u*∗, any errors in determining the latter will result in correlated variations in *C*_{D10n} and *ζ.* Where the inertial dissipation method has been used to determine *u*∗, these correlated variations have previously been incorrectly ascribed to an imbalance between production and dissipation of turbulent kinetic energy. By adding random variations to a dataset of calculated power spectral density values, the characteristics of this apparent imbalance have been determined. It is of similar magnitude to the imbalance term proposed by, for example, YT96 and DTWK97 (≈0.6*ζ*), and also shows the wind speed–dependent variations noted by YT96. Attempting to correct for this apparent imbalance will positively bias the calculated *C*_{D10n} values. However, compared to the scatter in most *C*_{D10n} datasets the effect is relatively small, about 10^{−4}.

For most wind speeds, provided the apparent imbalance is not corrected, random noise in determining *u*∗ will introduce scatter in the calculated *C*_{D10n} but no mean bias. However, in lighter winds under unstable conditions two effects can positively bias the calculated *C*_{D10n} values. Erroneously high *u*∗ will, because of the effect on the iterative calculation of *ζ,* bias *C*_{D10n} to a greater extent than erroneously low *u*∗ values. Also, for unstable low wind speed data with erroneously low *u*∗, there is a tendency for the iteration to not converge. Neglect of the nonconverged records positively biases the mean *C*_{D10n}. For typical datasets these effects appear at wind speeds below about 7 m s^{−1}. For relatively clean *u*∗ data this threshold wind speed will be lower (and vice versa). These low wind speed errors can be avoided by modifying the method of calculating *ζ.* Rather than using the experimentally determined *u*∗ value, a value calculated from a mean *C*_{D10n}–*U*_{10n} relationship is used. We have used a modified Smith (1980) formula, but the choice is not critical since it is only used to determine the stability corrections. Correcting the *C*_{D10n}–*U*_{10n} values published in Y98 for the erroneous inclusion of an imbalance term, and for the neglect of an anemometer response correction, resulted in a decrease of *C*_{D10n} by about 0.05 × 10^{−3}, which is probably not significant.

Understanding the cause of the apparent variation of *C*_{D10n} with stability increases our confidence in the inertial dissipation method. The potential importance of the assumed imbalance term was not that it changed the mean *C*_{D10n}–*U*_{10n} relationship significantly since much open ocean data is near neutral. Rather there was the implication that there was a stability dependence either of the balance between production and dissipation of TKE or of the sea surface roughness, or, more radically, that the similarity theory developed over land was not adequate over the ocean. Thus it appeared difficult to compare datasets taken under different stability conditions, for example, for offshore and onshore winds. Also, because the apparent stability effects varied with the quality of the data, different datasets appeared to show different stability dependence. By adopting the bulk stability estimate suggested here, these effects can be removed, significantly decreasing the scatter in the *C*_{D10n} estimates.

Finally we emphasize that most of the effects described here will occur in any dataset of experimentally determined *u*∗ estimates, not just those obtained from the inertial dissipation method. Where an increase of *C*_{D10n} has been observed as *U*_{10n} decreases below, say, 10 m s^{−1}, at least part of the increased *C*_{D10n} can be ascribed to random noise biasing the data. Indeed our own experience is that anomalously large *C*_{D10n} begin to be observed at higher wind speeds, say, the 6–10 m s^{−1} ranges, in data obtained from ships, whereas in datasets obtained from well exposed sensors on buoys the effects are only seen for winds below 5 m s^{−1} or less (YT97).

## Acknowledgments

This work was partially supported by the Joint Grant Scheme project, Coastal and Open Ocean Wind Stress. We thank Prof. K. L. Davidson of the Naval Postgraduate School Monterey for the use of the NPS data from HEXMAX.

## REFERENCES

## Footnotes

*Corresponding author address:* Dr. Peter K. Taylor, Southampton Oceanography Centre, European Way, Southampton SO14 3ZH, United Kingdom.