Variable and conflicting shear stress estimates inside a boundary layer with sediment transport

This paper presents a comparison between two methods for estimating shear stress in an atmospheric internal boundary layer over a beach surface under optimum conditions, using wind velocities measured synchronously at 13 heights over a 1.7 m vertical array using ultrasonic anemometry. The Reynolds decomposition technique determines at‐a‐point shear stresses at each measurement height, while the Law‐of‐the‐Wall yields a single boundary layer estimate based on fitting a logarithmic velocity profile through the array data.


Introduction
Shear stress (τ), or vertical momentum flux, is a crucial physical variable for understanding the dynamics of sediment transport by fluids (water or air) in turbulent boundary layer flows over the Earth's surface. It is traditionally deduced from the slope of the Prandtl-von-Kármán logarithmic velocity profile law (the 'Law-of-the-Wall').
Contemporary research into sediment transport remains grounded in the seminal works of Shields (1936) and Bagnold (1941) who developed predictive equations for sediment transport, typically relating sediment flux to the 1½ power of the time-averaged shear stress. Relating transport to time-averaged shear stress, however, often fails to accurately describe observed sediment flux in natural boundary layer environments, as it does not incorporate the effect of turbulent and unsteady flow conditions and related sediment flux response. Presently, there are no commonly used sediment transport equations that incorporate an explicit turbulence parameter.
In part, the predictive deficiencies may be attributed to methodological or technological limitations in measuring shear stress in natural boundary layer environments (Heathershaw and Simpson, 1978;Biron et al., 1998;Noss et al., 2010;Lee and Baas, 2012). Across the Earth sciences, however, improvements in sensor technology have provided methods for determining shear stress directly from turbulence statistics measured within the flow. For example, in hydrodynamic studies, acoustic doppler velocimetry (ADV) offers unobstructed three-dimensional flow vector measurements at high frequency in a small sampling volume (Voulgaris and Trowbridge, 1998) that enable shear stress estimates directly from flow turbulence, as a 'Reynolds' shear stress. Likewise, ultrasonic anemometry has provided a significant contribution to aeolian research by facilitating high frequency monitoring of wind velocities in three dimensions, and therefore the ability to derive Reynolds shear stress directly in atmospheric boundary layer flows. Furthermore, ultrasonic anemometers have become more robust and affordable, increasing their suitability for geomorphological applications (Walker, 2005).
The measurement of shear stress and its relation to sediment transport over a surface is typically framed within the structure of an internal boundary layer (IBL), the layer of fluid flow immediately above the ground where the flow dynamics are adjusted to the roughness and the no-slip condition at the contact surface. Terminology across meteorology, engineering, and Earth sciences varies considerably, with IBL, 'surface layer', 'logarithmic layer', 'inner region', and 'constant-stress layer' sometimes used interchangeably. We follow the definitions of Rao et al. (1974), where the IBL is composed of an inner equilibrium layer (IEL), typically the bottom 10% of the IBL, and a transition layer above. The logarithmic velocity profile is exhibited in a large part of the IBL (including the IEL), while a 'wake function' correction is needed in the upper parts of the IBL where the velocity profile is also partially controlled by 'outer region' scaling of the depth of the IBL. The IEL, however, is the region of the IBL that is considered the constant-stress layer (Garratt, 1990;Kaimal and Finnigan, 1994), where shear stress varies by less than 10% (Stull, 1988).
Under neutral convective stability (where thermally driven turbulence is negligible compared with mechanical shear turbulence) and uniform flow, the Law-of-the-Wall and Reynolds-derived shear stress estimates should theoretically be equivalent within the constant-stress layer (Tennekes and Lumley, 1972), nevertheless, empirical agreement between the two methods has been elusive. For example, in an aeolian boundary layer with no roughness elements Li et al. (2010) found that Reynolds-derived estimates were on average only~69% of the Law-of-the-Wall derived estimates; over complex terrain in wind tunnel simulations King et al. (2008) found that Law-of-the-Wall estimates overestimated Reynolds-derived estimates by an average 43%; and in a water flume Biron et al. (2004) found Reynolds-derived estimates were only~45% of those obtained from a logarithmic velocity profile. Ambiguous and conflicting estimates of shear stress are especially problematic for understanding sediment transport dynamics, as the 1½ exponent used in most predictive equations magnifies even small discrepancies. This paper aims to ascertain the relationship between Law-of-the-Wall (LOW) and Reynolds-derived (REY) shear stress estimates at multiple heights in the constant-stress layer of an atmospheric boundary layer using ultrasonic anemometry.

Background
Law-of-the-Wall shear stress Under conditions of neutral atmospheric stability and an aerodynamically rough surface, the vertical structure of horizontal wind speed in the lower levels of an internal boundary layer can be described by the Prandtl-von-Kármán logarithmic velocity profile law: where κ is the dimensionless von Kármán constant, typically 0.4, z is the measurement height (m) above the surface, z 0 is the aerodynamic roughness length (m) and U * the shear velocity (m s -1 ). The U * derived from the Prandtl-von-Kármán law may be estimated from the slope of the velocity profile when mean horizontal wind speed (U) at different elevations (z) is regressed against ln(z). A least-squares linear regression through these data points yields an estimation of: U * LOW = κ m, where m is the slope coefficient. The shear stress is then defined as: where ρ is the fluid density (kg m -3 ). The roughness length can be estimated from the same least-squares regression as: Over aerodynamically rough surfacessuch as those covered by dense vegetation elements or a well-developed saltation cloudthe physical ground surface may not accurately define the reference height datum for the Law-of-the-Wall, as the presence of roughness creates a so-called zero-plane displacement, or an (apparent) plane of momentum absorption lying at some small height above the ground surface. The zero plane displacement length (d) is incorporated in the Prandtlvon-Kármán logarithmic velocity profile law as an effective height adjustment: The displacement length depends on the roughness characteristics and must be determined from the wind profile, usually by numerical iteration over a fine-scale range of lengths to find the displacement that yields the best least-squares regression fit for the Law-of-the-Wall.

Reynolds derived shear stress
Shear stress may be calculated from the vertical turbulent momentum flux, using Reynolds decomposition to separate the wind vector into a time-averaged mean (denoted by an overbar) and a time-dependent fluctuating quantity (denoted by a prime) so that: U = Ū + U '. Three-dimensional ultrasonic anemometers, capable of resolving high frequency streamwise (u), spanwise (v) and vertical (w) wind velocity components, allow the turbulent momentum fluxes to be measured in the atmospheric boundary layer. If time-averaged velocities are aligned with local streamline coordinates (Stull, 1988), the streamwise Reynolds shear stress (one element of the Reynolds stress tensor) may be calculated as: τ ¼ Àρ u'w', where u'w' is the covariance between streamwise and vertical velocity components over a defined time period. Although this is the standard formula presented in the literature, to take into account the total shear acting along a horizontal plane the resultant horizontal Reynolds shear stress (τ_ REY ) should be calculated: where v'w' is the covariance between the spanwise and vertical velocity components.

Field experiment
A field experiment was conducted on the afternoon of 4 November 2011 at Tramore beach near Rosapenna, County Donegal, Ireland, on the edge of the Atlantic Ocean ( Figure 1). The beach is west-facing, approximately 4.5 km long and up to 150 m wide. The surface is un-vegetated with well sorted, quartz based, fine sand (0.19 mm mean diameter). The wind was from a south-south-westerly direction on an average bearing of 192°, staying inside a wind sector of 185°to 197°f or 90% of the time. An instrument array was installed on the beach surface on the high tide line 25 m from the toe of the foredune (N 55°09′50″, W 7°49′50″). The resulting fetch length under the mean wind direction was approximately 200 m. Following the recommendations of Jegede and Foken's (1999) field study of IBL structure the height of the IEL at the instrument array was estimated as: δ = 0.3(x) 0.5 , where x is the fetch length (200 m), yielding an IEL height of δ = 4.2 m. Three-dimensional wind velocities at 13 heights were measured synchronously in a close vertical array using horizontalarm ultrasonic anemometers (Gill model HS-50, velocity resolution: 0.01 m s -1 , accuracy: <± 1% RMS) operating at 50 Hz across a sonic path length of 0.134 m. The anemometers were aligned pointing into the mean wind direction and parallel with the beach surface, alternately mounted on adjoining masts, with the arms positioning the sensing volumes approximately 1 m upwind of the masts, establishing measurement heights at the centre of each sonic sampling volume of: 0.115, 0.230, 0.335, 0.460, 0.575, 0.680, 0.920, 1.020, 1.160, 1.265, 436 Z. S. LEE AND A. C. W. BAAS 1.375, 1.505, and 1.620 m ( Figure 2). These sonic anemometers and the accompanying data acquisition system have an established provenance of successful deployment in a number of other field studies (Lee and Baas, 2012;Delgado-Fernandez et al., 2013;Jackson et al., 2013;Smyth et al., 2014). Sediment transport was recorded simultaneously using 12 Safires of the same design outlined by Baas (2004). These sensors detect the impacts of saltating grains on an omnidirectional piezo-electric ring element 0.02 m in height, at an internal temporal resolution of 12.5 kHz, output and recorded at a frequency of 25 Hz. The Safires were positioned along a transverse array next to the anemometers, at 0.1 m intervals and with the sensitive ring 0.02 m above the sand surface.

Methodology
Data pre-processing From several multi-hour experimental runs, a 50 min sequence of data with periods of active sand transport was selected for further analysis. Quality control of wind data involved a series of pre-processing steps: (1) data from all sensors were synchronised within an accuracy of 0.016 s; (2) data were cleaned to replace occasional single-data-point anomalies (extreme, singular accelerations within the time-series) with the average of neighbouring samples; and (3) data were cleaned to remove spectral anomalies (extreme outliers to the frequency-size distribution of scalar speeds), replacing these with empty gaps in the time-series. The percentage of data removed amounted to less than 0.63%.
The 50 min sequence of data from each anemometer were rotated onto a standard right-handed x, y, z local coordinate system using a two-step streamline correction routine following Lee and Baas (2012). First, a yaw correction aligned the data according to the average wind flow vector in the x-y plane, so that, over the whole measurement run, v is zero. Second, a pitch correction rotated the data in the x-z plane to reduce, over the whole measurement run, w to zero. Resultant horizontal wind speeds at each height, , were calculated to produce estimates of τ_ LOW .
Safire data were synchronised to an accuracy of 0.032 s, and the signals normalised by dividing by the standard deviation of  the above-zero values in the time-series, following the procedure outlines in Baas (2008). The normalised signals were then averaged across the array to yield a single time-series of relative sand transport magnitude.

Data analysis
Shear stresses were determined at two temporal scales. First, for consistency with established literature, the total measurement run (a block of 50 min) was used as a single averaging period for the calculation of mean horizontal wind speeds Ū z , for τ_ LOWblock estimation and for the calculation of covariances for τ_ REYblock estimation.
Second, height-dependent moving-window averaging periods were used to generate time-series of both types of shear stress. The averaging period should relate directly to the eddy structure of the local fluid flow, which characteristically changes with height above the surface. An appropriate scaling parameter for turbulence analysis is the largest period (or largest eddy) associated with the local inertial sub-range, where kinetic energy is neither produced nor dissipated but is passed down to increasingly smaller scales (the eddy cascade) (Stull, 1988), the part of the spectral domain where power density is proportional to f -5/3 .
To define these 'bespoke', height-specific averaging periods, the inertial sub-range at each measurement height was identified with the second-order structure function of the local horizontal wind speed (Frisch, 1995;Martin et al., 2013) . A time-scale associated with the top of the inertial sub-range at each height, and thus the largest eddies, is deduced by identifying a break in slope in the structure function that corresponds with Δt 2/3 . For statistical robustness the averaging period should cover at least 12 of the largest eddies, and so the identified break-points were multiplied by 12 to yield averaging periods of 4. 91, 10.86, 14.09, 16.81, 18.73, 20.17, 22.76, 23.65, 24.75, 25.49, 26.21, 26.99 and 27.62 s, from bottom to top measurement height, respectively.
These periods were applied as moving-averages to produce smoothed time-series of horizontal wind speed at each elevation, which were used to derive a time-series of τ_ LOWwindow . Confidence intervals (95%) around the τ_ LOW estimates (on both time-scales) were calculated following Wilkinson (1984;see Namikas et al., 2003).
Time-series of τ_ REY for each measurement height were obtained by calculating covariances, over a moving window using the height-specific averaging periods, including a 'timelocal' streamline correction at these time-scales. 95% confidence intervals based on Heathershaw and Simpson (1978) were established using the standard error of the mean around u'w' and v'w': where N is the sample size.

Shear stress comparisons
On the time-scale of the whole measurement run (a 50 min block) τ estimated from the Law-of-the-Wall with a goodnessof-fit R 2 of 0.990, τ_ LOWblock = 0.061 N m -2 , is significantly higher than the 13 Reynolds-derived estimates for each measurement height, τ_ REYblock , as shown in Figure 3, most of which furthermore fall outside the confidence bounds of τ_ LOWblock . A particularly interesting result is that the τ_ REYblock estimates vary irregularly across the 'constant stress' boundary layer without a clear trend or pattern, and their variability greatly exceeds the confidence intervals around their individual estimates. The 95% confidence interval for τ_ LOWblock is 0.008 N m -2 , which represents 13.2% of its magnitude. The confidence intervals for τ_ REYblock are significantly smaller, trending from 0.0008 N m -2 (2.3% of the magnitude) at the lowest height to 0.002 N m -2 (4.6% of the magnitude) at 1.62 m.
The time-series of τ from both methods are compared by calculating overall run-means. As expected, for the τ_ LOWwindow time-series, the run-mean is identical to τ_ LOWblock , although the average confidence interval is larger (representing 21.8% of its magnitude). The run-means for the Reynolds-derived τ time-series, τ_ REYwindow , remain lower than the τ_ LOWwindow run-mean (except for height 1.505 m); however, they appear to show a slight trend of increasing with height approaching the magnitude of the τ_ LOWwindow run-mean at the top. At elevations above 0.575 m the τ_ REYwindow estimates, calculated using height-specific periods, are typically found to exceed τ_ REYblock , while below this elevation the difference is negligible. The confidence intervals, however, are significantly larger around τ_ REYwindow estimates, generally amounting to 27.3% of the shear velocity magnitude, averaged across all heights. Figure 4 shows the Law-of-the-Wall through the blockaveraged horizontal wind speeds. The R 2 of 0.990 demonstrates a linear fit that is excellent, particularly considering it involves 13 data points, confirming that this lower part of the internal boundary layer, i.e. the internal equilibrium layer (see above), is well developed and fully satisfies the Prandtl-von-Kármán logarithmic velocity profileas expected given the long uniform and flat fetch. The great and irregular variability in measured τ_ REY with height inside this IEL is, therefore, both significant and unexpected, since by definition this should be a flow zone of constant shear stress. Concerns about the possible effects of the masts and equipment potentially interfering with the local flow at the upwind measurement volumes of the sonic anemometer heads have been considered in detail. First, pairwise correlations of horizontal wind speed between neighbouring anemometers, at the original 50 Hz time resolution, nearly all greatly exceed 0.99 (Table I), suggesting a high degree of internal flow consistency that would be unlikely if there was a distorting effect of the downwind array structure, as this would involve flow divergence and significant deterioration of at least some, if not most, 50 Hz wind speed correlations. Second, wind tunnel studies of airflow through simulated porous shrub structures, which seem the most comparable with the sparse mast structure of our experimental set-up, suggest that any upwind projection of flow deflection or disturbance is limited to only a small distance away from the porous structure (Dong et al., 2008), and that in our array set-up the horizontal arms of the sonic anemometers are extending sufficiently upwind, roughly 1 m, to escape any such potential effects. Third, the exceedingly good fit of the Law-of-the-Wall, with an R 2 of 0.990, which is very high considering it involves 13 measurement heights, would be unlikely if the sensing volumes had been affected by flow obstruction. The only other field study with a comparable number of measurement heights is that of Namikas et al. (2003), who used eight cup-anemometers along a single mast (with little concern for flow interference) over a flat beach surface, and also reached goodness-of-fit levels on the order of 0.99.
In addition, the potential for instrument error or ultrasonic signal interference ('cross-talk') due to the close configuration of neighbouring anemometers can be excluded (Personal Communication: Gill Instruments). The significant differences between τ_ LOW and τ_ REY are also encountered in river channel flow studies where Law-of-the-Wall estimates are typically found to be higher than τ_ REY (Biron et al., 2004). Likewise, this trend was observed by Kim et al. (2000), working in the constant stress layer within a tidal boundary layer. Our results show that on the time-scale of the whole measurement run τ_ REYblock is equal to 68.0% of τ_ LOWblock , on average across all heights, and when based on time-series τ_ REYwindow equals 76.3% of τ_ LOWwindow . In comparison, Biron et al. (2004) found that Reynolds derived estimates were only 45.4% of those obtained from a logarithmic velocity profile. Kim et al. (2000) found that the estimates obtained using the Law-of-the-Wall were susceptible to the effects of boundary layer stratification due to suspended sediment loading (Smith and McLean, 1977) and variability in the depth of the constant stress layer. In our study, however, there is no evidence to suggest that the results presented here are a result of a poorly developed boundary layer, given the very high goodness-of-fit of the Law-of-the-Wall (with 13 measurement heights), as detailed above. Furthermore, since sediment transport in our field experiment occurred principally by saltation as opposed to suspension, and the saltation layer is limited to a height of 10 to 20 cm above the ground at best, we see no reason to suspect that some kind of boundary layer stratification is causing the discrepancies between our estimates. Nevertheless, for completeness we evaluate this possibility further below.
The inconsistency and irregularity in Reynolds-derived shear stress estimates within the constant stress layer are similar to those reported by Heathershaw and Simpson (1978), who present Reynolds stresses for two measurement heights, 0.5 m apart, in the bottom boundary layer of a macro-tidal current, with a coefficient of variation (CoV) of 13.7%. Results show a standard deviation over the 13 estimates of τ_ REYblock of 0.007 N m -2 , with a mean of 0.042 N m -2 , yielding a CoV of 16.0%. For run-mean τ_ REYwindow estimates, disregarding any possible trend, with a standard deviation of 0.009 N m -2 and a mean of 0.046 N m -2 , the CoV increases to 19.8%.
These unexpected findings are corroborated by analysis of a secondary dataset collected from a different 50-min sonic anemometry measurement run at the same Tramore field site on the same day (a dataset that was initially dismissed because it has no accompanying sand transport data). The same data quality assurance and analysis procedures were applied to yield shear stress estimates on the temporal scale of the whole measurement run (as a single block), shown in Figure 5.
The results corroborate the significant variability of Reynolds shear stress estimates with height, although they are more evenly distributed around the Law-of-the-Wall estimate and within the latter's confidence bounds. The latter are wider here than for the main measurement run because of a lower goodness-of-fit (R 2 = 0.943 through 11 points) associated with the Law-of-the-Wall estimated shear stress of τ_ LOWblock = 0.095 N m -2 Figure 4. Law-of-the-Wall fit through mean horizontal wind speeds at 13 heights, based on the whole measurement run as a single block, with grey-shaded 95% confidence intervals around the linear regression. This figure is available in colour online at wileyonlinelibrary. com/journal/espl (with 95% confidence intervals of 0.035 N m -2 ). The average of the Reynolds shear stress estimates over 11 measurement heights is 0.091 N m -2 , with a standard deviation of 0.017 N m -2 , yielding a CoV of 18.2%, comparable with the findings of the main run. These results further reveal that there is no consistency in the Reynolds shear stress variability between the two measurement runs at each height (in terms of over/underestimating relative to the average), suggesting that the Reynolds shear stress variations are related to the mechanics of the boundary layer flow rather than attributable to individual instruments.
Ambiguous estimates of shear stress are especially problematic for understanding aeolian sediment transport dynamics, as typically most predictive transport equations are based on the 1½ power of the time averaged shear stress (or the cubic power of U * ). For example, using Bagnold's (1941) predictive formula on τ_ LOWblock of the main run in our study yields a sand transport rate of 7.87 kg m -1 h -1 , whereas τ_ REYblock at z = 0.680 m above the surface predicts 3.38 kg m -1 h -1 (i.e. only 43%).
The inconsistency and irregularity in Reynolds derived shear stress estimates presents a significant problem for deciding the appropriate height for point-measurements of shear stress to predict sediment fluxes, particularly in studies where only one or two instruments are available. In fluvial studies, it is typically agreed that single point measurements should be taken outside of the roughness sublayer, but as close as possible to the bed surface (Babaeyan-Koopaei et al., 2002). Biron et al. (2004) suggest that single point measurements should be taken at a height equivalent to 0.1 of the flow depth, as this is the point where a maximal Reynolds stress is found in their work. Our results though show highest τ_ REY generally at the top of the profile, at elevations of roughly 1.5 m, or~35% of the IEL. Shear stress measurements at such heights above the surface are likely less relevant, however, for relating to aeolian sand transport dynamics that are largely restricted to within 10-20 cm above the bed, particularly in the context of high spatio-temporal variability (Baas and Sherman, 2006;Barchyn et al. 2014), and it is, therefore, probably preferable to measure closer to the surface rather than farther (see also next section below).

Time-series analysis
While the run-mean shear stresses are relatively low, significant temporal variability within the time-series of τ_ LOWwindow is observed, with estimates ranging from 0.006 N m -2 to 0.230 N m -2 . On average across all measurement heights, τ_ REYwindow ranges from 0.001 N m -2 to 0.154 N m -2 . In order to assess sand transport modelling impacts a saltation threshold shear stress (τ_ crit ) of 0.050 N m -2 is estimated using Bagnold's (1941) , where A is a constant (0.1 for fluid threshold), ρ s is the mineral grain density (2660 kg m -3 ), g is the gravitational acceleration (9.81 m s -2 ), and D is the mean grain diameter (0.19 mm). τ_ LOWwindow exceeds this threshold for 63.7% of the time, whereas τ_ REYwindow is above-threshold 34.6% of the time on average across all heights, ranging from 15.8% at the lowest elevation to 54.8% at 1.505 m. This presents a difficult comparison with the level of sand transport recorded by the Safires at 47% of the total measurement run, although the transport threshold may be better conceptualised as a stochastic variable, which would yield variability in the exceedance fractions that might match better with the recorded transport activity Throughout the main measurement run, the goodness-of-fit for τ_ LOWwindow is consistently above 0.900 (95.7% of the time) and frequently above 0.950 (81.2%). Considering the relatively short averaging windows used to smooth the velocity data and derive τ_ LOW , and given that the linear fit is regressed through 13 heights, the high percentages of time-series results exceeding 0.900 and 0.950 R 2 is encouraging.
Time-series of τ_ LOWwindow and τ_ REYwindow correspond poorly with each other throughout the whole measurement run, as exemplified by Figure 6(a) (which shows the Law-ofthe-Wall time-series compared with τ_ REYwindow at a height of 1.020 m, as representative of typical field studies). The two methods often yield both conflicting estimates as well as opposite trends, e.g. between 420 and 440 s. On a time-variable basis, the instantaneous disparities between the methods are very wide-ranging, with height-averaged τ_ REYwindow ranging between 2.8% and 1258.7% of τ_ LOWwindow . Furthermore, this poor correspondence is confirmed by low correlation coefficients between the time-series of τ_ LOWwindow and τ_ REYwindow at each measurement height. All correlations are found to be lower than 0.41, and no clear sequence is exhibited with elevation.
The disparities between individual τ_ REYwindow estimates at different elevations on a time-variable basis far exceed those reported above for the whole measurement run. The CoV averages to 39.3% over the time-series, but ranges between 9.6% and 244.8%.
The goodness-of-fit of the Law-of-the-Wall (Equation (1)) may be improved by incorporating a displacement length, d, to offset height z (Equation (2)) when the momentum absorption plane (the zero-plane) is not equivalent to the datum employed in the field (the physical sand surface). This adjustment may be required to correct the elevations of the horizontal wind speed measurements to a single time-averaged displacement of the momentum absorption plane, or in the case of fluctuating sediment transport the adjustment may require a time-variable zero-plane displacement. The former was explored by testing for improvements in R 2 when deriving τ_ LOWblock estimates using a range of displacement lengths between -0.115 and +0.115 m (the height of the lowest anemometer) at intervals of 0.001 m. The highest R 2 was found with a displacement length of zero, confirming that for the measurement run as a whole there was no consistent zero-plane displacement. The same systematic approach was used to evaluate if a time-variable displacement length should be employed in the time-series of τ_ LOWwindow , to potentially account for saltation activity (Figure 6(b) and (c)). Our analysis shows, however, only a minimal improvement in the goodness-of-fit, raising only an additional 0.5% of the data above an R 2 of 0.900 and an additional 4.5% above the 0.950 threshold. More importantly, as shown in Figure 6, the variable displacement does not appear to relate temporally to levels of sand transport activity and furthermore the length magnitudes seem unrealistic, reaching heights of up to 0.1 m above the bed.
The vertical variability in Reynolds-derived shear stress found in our study presents a fundamental problem for relating wind forcing to sand transport, and so as part of time-series analysis we explored the temporal correlations between sand transport and (τ_ REYwindow ) 1.5 for the 13 different heights. The time-series of shear stress and sand transport were blockaveraged to a 5 Hz temporal resolution and analysed for cross-covariance as a function of lag time, to yield the maximum cross-covariance at optimum lag. The optimum lag increases for shear stress measured further away from the surface (not further discussed here), but more importantly the associated covariance increases for shear stress measurements closer to the bed, as shown in Figure 7. The results display a small, but statistically significant, vertical trend with a correlation coefficient of -0.900, suggesting that for sand transport studies τ_ REY should be measured as close to the surface as conditions permit.
Exploring potential sediment transport effects While sand transport does not appear to justify any timevariable zero-plane displacement, the presence of active saltation may still somehow be a potential cause for the disagreement between τ estimates. First, saltating grains have the potential to interfere with the lower anemometers directly by interrupting the sonic pathway, or impacting on the ultrasonic transducers. Second, it has been argued that sand transport may indirectly affect estimates of τ_ LOW due to sediment loading changing the fluid properties of the flow, necessitating a 'variable' von Kármán constant for Law-of-the-Wall calculations (Smith and McLean, 1977;Li et al., 2010).
The potential effect of saltation is evaluated by comparing estimates of τ derived during episodes of active sand transport against periods without activity. To allow comparisons between τ_ LOW and τ_ REY across all measurement heights, the pre-

CONFLICTING SHEAR STRESS ESTIMATES INSIDE A NATURAL BOUNDARY LAYER
processed synchronous data of both wind and sand transport were first block-averaged using the relevant height-specific averaging period. Each block was then ranked according to the average relative magnitude of transport recorded by the Safires, and 'transport' episodes were selected as those highestmagnitude blocks that together comprise a cumulative 90% of the total transport. These were matched against an equal number of randomly sampled inactive 'clean air' blocks.
Using these two suites of episodes, Figure 8 demonstrates that irrespective of the presence or absence of sand transport τ_ REY exhibits a similar vertical variability to the total run-mean estimates of τ_ REYblock (cf. Figure 3). The clean air episodes exhibit a standard deviation of 0.007 N m -2 in τ_ REY over 13 heights, with a mean of 0.034 N m -2 , yielding a CoV of 18.9%. Similarly, τ_ REY during transport episodes shows a standard deviation of 0.014 N m -2 , a mean of 0.063 N m -2 and a CoV of 22.5%.These results are comparable with the run-mean τ_ REYblock CoV of 16.0%. This suggests, therefore, that the variable and conflicting estimates of Reynolds shear stress are not associated specifically with sand transport activity.
Scatter plots for each measurement height comparing τ_ LOW and τ_ REY differentiated according to the two types of episode are shown in Figure 9. The association between clean air and transport episodes depends on measurement height: close to the surface there is considerable overlap of the two data groups, whereas at higher elevations they diverge, with transport episodes occurring under clearly stronger shear stresses. Linear regressions through both types of data groups, with forcing through the origin, show that close to the surface regression slopes are typically steeper for both types of episodes, indicating that nearer the surface τ_ REY appears to increasingly underestimate compared with τ_ LOW , as was already apparent in Figure 3. Note however, that nearer the surface, the regression slopes for clean air exceed those for transport episodes,  particularly at the three lowest elevations. For both types of episodes, the slopes start to approach unity at elevations of 0.46 m and above.
With the exception of the lowest two measurement heights, the similarity between regression slopes for both transport and clean air episodes, particularly in the context of their 95% confidence intervals, suggests that there is no sand transport impact on the relationship between τ_ REY and τ_ LOW , and a transportdependent 'apparent' (or 'variable') von Kármán constant is not suitable for relating the two. This contrasts with hydrodynamic studies which observe a stratified boundary layer induced by sediment loading, causing a divergence between at-a-point vertical turbulent momentum flux and τ determined from a velocity profile (Smith and McLean, 1977;Adams and Weatherly, 1981;Kim et al., 2000). Li et al. (2010) postulate that this effect is responsible for incompatible τ_ LOW and τ_ REY estimates at 1 m above the bed, during periods of active aeolian sand transport in their field experiments. Density stratification and the associated 'apparent' von Kármán constant are unlikely over the larger part of the vertical profile in our study, however, as our results show, only the lowest two measurement heights demonstrate any significant divergence in regression slopes between clean air and transport episodes. The observed differences between slope estimates at the lower elevations are, furthermore, the exact opposite to the postulated expectations of sediment-laden boundary layer stratification, as our results nominally suggest that during active transport episodes the von Kármán constant would require less adjustment to match τ_ LOW to τ_ REY than under clean air conditions.

Evaluation
In boundary layer flows τ_ REY and τ_ LOW are typically assumed to be equivalent within the constant stress layer, but the observations presented here indicate that they are fundamentally conflicting on both a run-averaged and a time-series basis. Our results furthermore show that τ_ REY varies irregularly across the 'constant' stress layer, with a 15-20% variability that greatly exceeds individual confidence intervals.
Data limitations or errors are unlikely sources for the discrepancies observed here: quality was assured through a careful experimental design and data selection, rigorous data cleaning, and robust pre-processing and analysis. Furthermore, spectral analysis of the 50 Hz data at 0.115 m (i.e. the lowest elevation) indicates that the Kolmogorov micro-scale is reached before the Nyquist frequency, suggesting minimal exclusion of high frequency turbulent energy from the data ( van Boxel et al., 2004).
It is acknowledged that sediment transport has the potential to cause disagreement between shear stress estimates directly, by interrupting the sonic pathway, or impacting the ultrasonic transducers. While the discrepancy between τ_ LOWblock and the τ_ REYblock at 0.115 m is significant (56.4% of the τ_ LOWblock magnitude), the difference is comparable with estimates derived at 0.68 m (56.9% of τ_ LOWblock ) or 1.16 m (59.2% of τ_ LOWblock ). The discrepancies are, therefore, unlikely due to mechanical interference of saltating grains. More importantly, our results indicate that there is no justification for an 'apparent' or variable von Kármán constant and we find that sediment transport is not likely responsible for the conflicting estimates of τ_ LOW and τ_ REY .
The discrepancies are not simply within normal error bounds, and the τ estimates are based on robust statistics. For example, we record an R 2 of 0.990 for the fit of the Law-ofthe-Wall (through 13 points), an R 2 comparable with many field and wind tunnel studies using far fewer measurement heights (Butterfield, 1999;Li, et al., 2010). Furthermore, our confidence bounds are equal to 4.6% of the shear stress, smaller than those reported by Namikas et al. (2003), for example, which range from 7.1% to 12.0%.
In summary, the observed discrepancies, both between the two methods, as well as between heights in terms of Reynolds derived stress, are significant and unlikely due to data or instrument error, nor due to sand transport effects. The discrepancies not only present significant problems for studying sediment transport dynamics in natural boundary layers, but our findings also challenge the notion that shear stress within the theoretical framework of the Law-of-the-Wall is equivalent to the shear stress derived from the local Reynolds stress tensor, as well as the presence of a 'constant' stress layer.
Our results show that a typical aeolian internal boundary layer flow as is commonly measured in field studies may display a near-perfect logarithmic velocity profile and yet at-apoint Reynolds shear stresses can vary considerably with height and deviate significantly from the Law-of-the-Wall estimate. The internal variability in τ_ REY of 19.8% for window-averaged estimates and 16.0% for block-averaged estimates is far in excess of the definition of a constant stress layer, where variability in fluxes by elevation should remain below 10% of their magnitude (Stull, 1988).
The results demonstrate that τ_ LOW is strictly a bulk property associated with the entirety of the flow, while τ_ REY corresponds to a single point within that flow, and critically the two are not equivalent on relevant time-scales. One potential cause for the discrepancies, we hypothesise, is that the presence of coherent flow structures in the turbulent flow (cf. Bauer et al., 2013) creates highly localised and persistent variability in the stress tensor at different points in the flow (hence affecting the Reynolds shear stress), but they have less impact on the scalar wind speed used in the Law-of-the-Wall framework. Variability in streamwise vorticity, for example, will affect the stress tensor but not the scalar wind speed. The potential influence of bursting events on flux profiles (Chen, 1990) may explain shear stress variability on a time-series basis, on the temporal scale of the events themselves, but would not explain the longer-term discrepancies and variability on the time-scale of the whole measurement run found in this study.
While we do not mean to challenge the validity of the constant stress layer concept (fundamentally a conservation of momentum principle), our findings do question the typical expectation of the presence of a constant stress layer under the kind of field conditions and site characteristics that we would normally consider to be ideal and representative for (aeolian) sediment transport studies (such as a long fetch distance, unidirectional winds, a flat and uniform sand surface, and a very high R 2 goodness-of-fit for the Law-of-the-Wall). Consequently, it is difficult to differentiate between 'simple' boundary layer conditions where Biron et al. (2004) recommend using Reynolds derived shear stress, or more complex flow conditions where Reynolds derived methods are found to be less reliable (Kim et al., 2000). The use of streamline correction routines (van Boxel et al., 2004;Lee and Baas, 2012) and the calculation of the resultant horizontal shear stress, however, should improve the use of Reynolds derived techniques irrespective of the complexity of the surface topography.
Further research comparing the two methods with a simultaneous direct measure of surface shear stress, using a drag plate in a full-scale field experiment, may yield firmer guidance on measurement strategy and application. At present we recommend using τ_ REY , albeit with careful consideration (streamline correction, resultant horizontal). Our preliminary analysis of correlations between τ_ REY at different heights and sand transport suggests that for studies concerned with evaluating sediment transport dynamics τ_ REY should be measured as close 443 CONFLICTING SHEAR STRESS ESTIMATES INSIDE A NATURAL BOUNDARY LAYER to the surface as possible, while avoiding any potential mechanical interference of saltating grains.
It may be fruitful to evaluate potential relationships between turbulent kinetic energy (TKE) and sediment transport to facilitate a move away from variable and conflicting estimates of τ. This suggestion is in agreement with recent studies of airflow over more complex surface topographies, for example: Smyth et al. (2014) concerned with flow dynamics and associated sediment transport patterns within a coastal trough blow out, and Chapman et al. (2013) who evaluated Reynolds stress and associated sediment transport over a coastal foredune. Methods to derive τ using a modified TKE are used in fluvial geomorphology (Biron et al., 2004;Pope et al., 2006;Noss et al., 2010) based on a linear relationship between TKE and τ (Soulsby and Dyer, 1981;Kim et al., 2000). Research comparing TKE with τ_ LOW and τ_ REY has typically found that the TKE method is more suitable in complex flow conditions (Biron et al., 2004;Bagherimiyab and Lemmin, 2013). While it is proposed that more specific investigation into TKE and sediment transport dynamics is required, the existing research on TKE derived shear stress estimation is encouraging.