1 Introduction

Understanding flow around urban environments is becoming increasingly important as cities and their populations grow in size (Department of Economic and Social Affairs 2019). Although surface energy balance models have recently improved, accurate models for aerodynamic parameters are still poor, particularly for non-conventional roughness geometries—e.g., tall canopies with heterogeneous height, with urban flows that remain poorly described by both theoretical and empirical models (Kanda et al. 2013). The surface layer within the atmospheric boundary layer (ABL) is customarily split into two regions. First, the roughness sublayer (RSL) which acts from the wall surface up to a certain point above the roughness elements. Customarily, the RSL is the region where the flow still encounters the effect of the individual roughness elements (Kent et al. 2017). Secondly, the inertial sublayer (ISL), which covers a region above the RSL (Raupach et al. 1991). The ISL begins at the top of the RSL and is characterized as a region of constant momentum flux. Considerable debate has taken place over defining the boundaries of both the RSL and ISL. The height of the RSL has been commonly quoted between 2–5 times the average height (h) of the roughness surface (Raupach et al. 1991), and more recently as low as 1.1–1.2h (Florens et al. 2013). The ISL has been subject to even more debate. Traditional theory (Stull 1988) suggests that, in a sufficiently developed boundary layer, an ISL should form. In contrast, Jiménez (2004) postulated that for \(\delta /h<80\) (where \(\delta \) is the boundary-layer height), the RSL increases in depth and effectively replaces the ISL. The Jiménez (2004) idea has been supported by several studies (Rotach 1999; Cheng and Castro 2002; Cheng et al. 2007; Hagishima et al. 2009). However, Leonardi and Castro (2010) and Cheng et al. (2007) argue that the relative boundary-layer depth quoted by Jiménez (2004) may not be accurate. Here, we also refer to a constant-stress layer (CSL) in Sect. 3.2.4, which is to be considered a relaxed version of the ISL region of Jiménez (2004) that does not necessarily require either the existence of a logarithmic velocity profile or the universality of the mean turbulence statistics.

Typically, in the standard characterization of the effect of wall roughness, a surface (at least close to the wall) is characterized by the aerodynamic parameters and the friction velocity. It is of increasing importance to be able to determine the zero-plane displacement (d) and the roughness length (\(z_{0}\)) to effectively predict and calculate the flow in and above the growing urban environments (as discussed by Stull (1988) in Sect. 9.7). In fully-rough conditions, these parameters are usually obtained, following Cheng and Castro (2002), by (i) calculating the friction velocity (\(u_{*}\)) and then (ii) applying a logarithmic law fitting procedure in the constant-flux layer

$$\begin{aligned} U=\frac{u_*}{\kappa }\ln \left( \frac{z-d}{z_0}\right) . \end{aligned}$$
(1)

Here, U is the wind speed, \(\kappa \) is the von Kármán constant, and z is the vertical coordinate. The friction velocity can be determined in two ways. Directly, by calculating the drag generated by the rough wall, for example using fully instrumented elements with static pressure ports or by a floating-element force balance (Cheng et al. 2007; Hagishima et al. 2009; Zaki et al. 2011). These methods are both acceptable in the fully-rough regime, where the viscous drag is almost negligible in comparison with form drag (Leonardi and Castro 2010). The second approach is indirect, where the friction velocity is evaluated from the Reynolds shear stress in the ISL (Cheng and Castro 2002; Cheng et al. 2007), so that \(u_{*}=\sqrt{\tau _{0}~\rho ^{-1}} \approx \sqrt{- \overline{u'w'}}\), where \(\tau _{0}\) is the wall shear stress, \({\rho }\) is the density of the medium, and \(u'\) and \(w'\) are the fluctuating streamwise and vertical velocity components. The roughness length \(z_{0}\) is the height at which the wind speed reaches zero, and is used as a scale for the roughness of the surface (Stull 1988). The zero-plane displacement d is understood as a correction factor to the logarithmic profile (Cheng and Castro 2002; Kanda et al. 2013), while physical reasoning tends to describe it as the central position at which drag occurs on a rough wall (Jackson 1981). The aerodynamic parameters vary due to surface properties and, therefore, it should be possible to examine their impact on the flow by systematically varying the surface characteristics. Many experiments have taken place over uniform-height cube arrays with different values of the parameters h, \(\lambda _{p}\), and \(\lambda _{f}\) (see Grimmond and Oke (1999) for the packing densities \(\lambda _{p}\) and \(\lambda _{f}\) definitions) to examine how roughness affects the aerodynamic parameters (Cheng and Castro 2002; Cheng et al. 2007; Jackson 1981; Grimmond and Oke 1999; Florens et al. 2013; Macdonald et al. 1998; Sharma and García-Mayoral 2020). However, the discrepancy in these studies, with additional evidence from some recent numerical and experimental work, have demonstrated the inaccuracy of using just these variables to examine the sensitivity of the aerodynamic parameters (Hagishima et al. 2009; Zaki et al. 2011; Kanda et al. 2013; Kent et al. 2017; Nakayama et al. 2011; Xie et al. 2008). Previous work has also highlighted how \(\lambda _{p}\) or \(\lambda _{f}\)—in isolation—are insufficient to characterize non-cubical roughness (Carpentieri and Robins 2015), and has advocated for the need to decouple the two solidity ratios (Placidi and Ganapathisubramani 2015, 2017); however, this is outside the scope of this work. Kanda et al. (2013) pointed out the importance of two additional parameters when describing a canonical regular surface, the standard deviation in height variations (\(\sigma _{h}\)), and maximum element height (\(h_{max}\)). Others have attempted to determine the aerodynamic parameters for various environments by deriving semi-empirical relationships based on roughness geometry (Kanda et al. 2013; Macdonald et al. 1998; Kent et al. 2017). These have the advantage of being able to quickly determine the parameters without the need for field observations, wind-tunnel tests, or computational experiments, allowing for the fast prediction of flow in urban environments. Here, we further explore how the aerodynamic parameters behave in extremely rough surfaces with heterogeneous heights in comparison with uniform-height cases at a matching average height. The block arrays are based on simplified attributes of super-tall grid cities, which feature a large standard deviation of element heights with large aspect ratios. This represents a combination of variables, which to the authors’ knowledge has not yet been explored. Below, the experimental set-up is discussed in Sect. 2. Sections 3.1.13.1.4 and 3.2.13.2.4 examine the depths of the boundary layer and surface layers. The aerodynamic parameters are then presented in Sects. 3.1.5 and 3.2.5. Finally, conclusions are drawn in Sect. 4.

2 Experimental Facility and Details

2.1 Experimental Facility

Experiments were conducted in the Aero tunnel within the EnFlo laboratory at the University of Surrey. This is a closed-circuit wind tunnel with a maximum flow speed of \(40~\hbox {m s}^{-1}\). The freestream flow speed, measured by a Pitot tube located upstream of the model was set to \(10~\hbox {m s}^{-1}\) for all cases presented here. The tunnel’s test section is 9 m long, 1.27 m wide, and 1.06 m tall. The streamwise, spanwise, and vertical directions are identified with the x, y, and z axes, respectively. The z-axis is set from the top of the baseboard of the model (i.e., the actual wall), while y = 0 is set in the centre of the test section. The position x = 0 is considered as the beginning of the tunnel test section. Time- and spanwise-averaged mean and fluctuating velocity components are denoted as (UVW), and (\(u',v',w'\)), respectively. In the space between the beginning of the test section and the model, a 1-m-long ramp rises from the floor of the tunnel to the average canopy height (\(h_{avg} = 80~\hbox {mm}\)). The ramp creates a smooth transition between the wind-tunnel test-section inlet and the roughness surface, thus minimizing the flow disruption at the beginning of the roughness fetch (Cheng and Castro 2002) and allowing an equilibrium boundary layer to form. The most upstream measurement station is at \(x = 3600~\hbox {mm}\), where the flow is already fully developed.

2.2 Rough-Wall Models

Four surface roughnesses, all representing idealized tall and super-tall urban environments, were used here. Two of the surfaces have elements of uniform height, while two have elements with varied height. The individual roughness elements are sharp-edged cuboids with an average height of 80 mm. Urban buildings are classified by their aspect ratio, \(AR=h/w\) (where h and w are the height and width of the building, respectively). Buildings with \(3<AR<8\) fall in the tall regime (generally \(100~\hbox {m}<h< 300~\hbox {m}\)), while buildings with \(AR>8\) (and \(h> 300~\hbox {m}\)) are referred to as super-tall (Council on Tall Buildings and Urban Habitat 2019). Based on this criterion, the surface morphologies examined here are classified as tall or super-tall. Zero-pressure-gradient conditions were used here, as the acceleration parameter, calculated for both uniform- and varied-element height cases, was \(4.85 \times \ 10^{-8}\) and \(9.48 \times \ 10^{-8}\), respectively. The surfaces investigated are further described in the following sections.

2.2.1 Homogeneous-Height Model

To create an idealized urban environment in the case of a uniform-height canopy, two urban features were studied: the packing density, \(\lambda _p\), and the element aspect ratio. Large urban areas with high-density buildings of 14 dense cities were used to calculate a characteristic urban packing density, as reported in Table 1. The Google Maps facility was used to measure the area plots and building base sizes. A range between \(\lambda _{p} = 0.33\) to 0.57 was determined, which is in line with values cited for real cities in Kent et al. (2017). The packing density of \(\lambda _{p} =0.44\) was within this range, and it was selected to describe densely packed cities. The Council on Tall Buildings and Urban Habitat (2015) calculates a mean height of all buildings in Manhattan (New York) above 100 m to be 145.7 m. Using this average height and the bases of buildings from Google (2021) gives an aspect ratio of approximately 3.4. Guided by this criterion, for the effects of walls and fetch, a scaled model with elements of size \(80~\hbox {mm} \times 20~\hbox {mm} \times 20~\hbox {mm}\) was designed. The roughness was mounted on five base plates. Each base plate when rotated allowed the roughness pattern to vary from aligned to 50% staggered, as in Fig. 2. These surfaces are referred to as uniform-height aligned (UHA) and uniform-height staggered (UHS). In the homogeneous-height roughness cases, a total of 5775 elements were used.

Table 1 Packing densities of highly dense areas within large urban areas
Fig. 1
figure 1

A single board of a homogeneous and b heterogeneous elements

Fig. 2
figure 2

Aligned and staggered orientations of \(3 \times 3\) module for homogeneous and heterogeneous models (all measurements are in mm). Indicated in the centre of the grey blocks are the heights of the elements in a heterogeneous module. Indicated with blue dots and black lines is an example of where the vertical profiles were measured within a repeating unit

2.2.2 Heterogeneous-Height Model

A varied-height model was also designed that differed from the uniform-height model in only one aspect: the standard deviation of the element height (\(\sigma _{h}\)). The two different roughness configurations are shown in Fig. 1, where the aligned homogeneous and heterogeneous surfaces are depicted in (a) and (b), respectively. The use of a large \(h_{avg}\) value provides the ability to introduce large \(\sigma _{h}\). The standard deviation of the elements was modelled on the districts of Mong Kok (Hong Kong) and Midtown Manhattan. The geometric properties of the elements were derived from the NYC Department of City Planning (2018) and Hong Kong Map Service (2019), along with the buildings’ design guide by Bull (2012), when the information was lacking in the database. The real \(\sigma _{h}\) of cities is 17.3 m and was scaled down to \(\sigma _{h} = 0.049~\hbox {m}\). As in the homogeneous model, the average height is 0.08 m, hereafter denoted as \(h_{avg}\). Other than \(\sigma _{h} = 0.049~\hbox {m}\), element heights were selected by matching the distribution of building heights in the dataset, resulting in an increased number of short elements, with only a few elements exceeding the value of \(h_{avg}\), where \(h_{max}=0.2~\hbox {m}\). The rough surfaces were constructed by assembling elements into modules. Herein a module is \(3 \times 3\) elements consisting of five different heights, randomly placed (Fig. 2). In the heterogeneous model, a repeating module is needed to achieve statistically relevant statistics, as described by Cheng and Castro (2002). Each module contains elements with heights 30, 30, 60, 60, 60, 80, 80, 110, and 200 mm, distributed as in Fig. 3. The purpose of this randomization is to avoid creating preferential corridors for the flow between the super-tall elements facing the wind direction and to create a more realistic city layout. The tallest elements in the varied-height surface are within the super-tall regime, with \(AR=10\) (Jianlong et al. 2014). The varied-height surfaces have a total of 5445 elements. Similar to the homogeneous canopy, each base plate when rotated allowed for the roughness pattern to be modified from aligned to 50% staggered. Due to additional geometrical constrains, the total number of elements is slightly different for the heterogeneous and homogeneous canopies. However, this is of little consequence as the downstream flow is spanwise self-similar over the different repeating units and fully developed in a central section of the wind tunnel (i.e., at the measurement stations). The surfaces are referred to as varied-height aligned (VHA), and varied-height staggered (VHS). A summary of the surfaces can be found in Table 2.

Fig. 3
figure 3

Probability density function of the distribution of element heights per module in the heterogeneous model

Table 2 Characteristics of the different roughness surfaces considered

2.3 Instrumentation

A two-component laser Doppler anemometer (LDA, FiberFlow, Dantec Dynamics, Denmark) was employed to measure two velocity components of the flow simultaneously. The laser beams were converged by a 300 mm focal length lens, which facilitated measurements in between elements. In the varied-height model, two elements of 200 mm height are located adjacent to each other in some configurations. In these situations, no measurements can be acquired below \(z = 40~\hbox {mm}\) due to laser obstruction by the elements. The model was sprayed with black matt paint to minimize light reflections. The LDA probe was mounted onto a traverse system which could move the measurement volume in three-dimensions within the tunnel with sub-millimetre accuracy. An elliptical mirror was attached to the LDA probe to rotate one of the laser-beam couples by \(90^\circ \), effectively changing the component of velocity measured by the LDA.

Additionally, a fully-instrumented pressure-tapped element was used in the case of the uniform-height surfaces to measure the differential force across an element, and thus its drag. To allow for these measurements, one element was removed and replaced with an identical three-dimensional-printed plastic element, fitted with a total of 25 static pressure ports on one of its faces. The element could then be rotated to allow for the differential pressure to be measured. The drag force was determined by integrating the pressure distribution over the front and back of the pressure-tapped element, allowing measurement of the shear stress due to pressure \(u_{*}(p)\) as in Cheng and Castro (2002). The value of \(u_{*}(p)\) was only determined for the uniform-height surface given the impossibility of replicating all random permutations of the different height elements for the heterogeneous roughness case.

2.4 Experimental Details

An LDA was used to measure velocity and Reynolds stress profiles within the ISL and RSL above the elements, but also within the canopy. For each surface, different types of LDA profiles were acquired as shown in Fig. 4. For both homogeneous and heterogeneous models, x-profiles were collected to study the development of the boundary layer. In the homogeneous model, nine vertical profiles were taken over a repeating element (see dots in Fig. 2) to evaluate the depths of the ISL and RSL. Over the varied-height surfaces, 81 vertical profiles were taken over one module to study the statistical convergence of the Reynolds shear stress profiles within the RSL (i.e., \(3\times 3\)-module in Fig. 4). Finally, over the heterogeneous model, 18 vertical profiles across the width of the tunnel were taken to study the CSL depth (further details on y-planes are contained in Fig. 2).

Fig. 4
figure 4

Schematic indicating the location of the LDA measurements. All units are in mm

3 Results and Discussion

This section presents results for the cases of uniform height in Sect. 3.1 followed by results for the varied-height cases in Sect. 3.2.

3.1 Effect of Uniform Element Height

3.1.1 Depth of the Boundary Layer

The depth of the boundary layer \(\delta \) is commonly estimated for the height where \(U = 0.99\) \(U_{ref}\). The freestream flow speed (\(U_{ref}\)) used in all plots is the flow speed measured at the wind-tunnel inlet, ahead of the rough wall models. Since the roughness length \(z_{0}\) is a strong function of roughness fetch in a developing flow (Cheng et al. 2007), the boundary layer must be in equilibrium with the surface below it and be fully-rough before representative measurements can be collected. It was determined that this occurs at a fetch of \(x = 4000~\hbox {mm}\) (\(\approx 15\delta \)) for both varied- and uniform-height models. The results past this location are shown in Fig. 5, where \(\delta \) is shown to be up to the height of 3.25h and 3h in the UHS and UHA cases, respectively. The increase in boundary-layer thickness in the UHS case is likely due to increased frontal blockage of the staggered array. When the elements are staggered, a ‘skimming’ flow regime (Grimmond and Oke 1999) is less likely to occur, as the streamwise distances between the elements increase, allowing for the development of the ‘wake interference’ flow regime, with associated enhanced turbulent structures and an increased boundary-layer thickness.

Fig. 5
figure 5

Boundary-layer depth shown with vertical profiles of streamwise velocity component over UHA and UHS models

The current data (\(h_{avg} = 80~\hbox {mm}\) and \(\lambda _{p} = 0.44\)) are compared with the previous literature in Table 3. The taller elements in this experiment occupy a much larger fraction of the boundary layer (\(\delta /h_{max}=3.25\)) than those of Cheng and Castro (2002) and Cheng et al. (2007), where \(\delta /h_{max}=12\) and \(\delta /h_{max}=7\), respectively (i.e., the surfaces described here have a higher relative roughness height). At the same freestream flow speed, the elements of \(h = 20~\hbox {mm}\) of Cheng and Castro (2002) and Cheng et al. (2007) show a much larger value of \(\delta /h_{max}\) than found here, indicating that—not surprisingly—the height of the elements can influence the boundary-layer depth. For the sake of brevity, the reader is referred to Andrieux (2017) and Thorpe (2018) for a more in-depth discussion of the boundary-layer variation across stations in both the streamwise and spanwise directions due to surface heterogeneity.

Table 3 Surface characteristics of various roughness configurations

3.1.2 Mean and Fluctuating Velocity Profiles

Mean streamwise velocity component profiles for flow over uniform-height cases are shown in Fig. 6a. These profiles have been spatially averaged across one repeating unit to offer a fair representation of the surfaces under investigation and are taken in the region of fully-developed flow. The units are non-dimensionalized with the skin friction velocity \(u_{*}\) and the roughness length \(z_0\). For both cases in Fig. 6a, a linear region within the profiles is evident; this offers strong support for the existence of a well defined ISL. Within this region, the collapse of the rough cases onto the smooth theoretical dashed line is an indication of the validity of the assumptions underpinning the methodology used to estimate the roughness parameters. This is further discussed in Sect. 3.1.5. The logarithmic region seems to extend well into what would be expected to be the roughness sublayer, as suggested in Cheng and Castro (2002). It must be stressed that data from different streamwise stations were considered for this analysis. The characteristics of the linear region visible in Fig. 6a, and the roughness parameters describing it, were found to be fully converged and independent from the streamwise location (i.e., from the boundary-layer depth) once the flow developed over approximately 15\(\delta \). This offers further proof that \(15\delta \) is the minimum required roughness fetch for the flow to be in equilibrium with the underlying surfaces considered herein. The parameters reported in Table 3 refer to these conditions, where the flow is fully developed. Next the streamwise turbulent fluctuations are discussed in Fig. 6b in the form of the diagnostic plot (Alfredsson and Örlü 2010), which removes the need to define the roughness parameters and the friction scaling. The value of \(U_{e}\) used to normalize the mean streamwise velocity component in Figs. 6 and 10 is the local velocity component at the edge of the boundary layer. The benefit of the diagnostic plot approach here is twofold: first, it allows one to assess the universality and the self-similar character of the turbulent statistics; second, it provides an independent assessment on whether the flow is in fully-rough conditions. For these reasons, the smooth (Alfredsson et al. 2011) and the fully-rough (Castro et al. 2013) asymptotes are reported in Fig. 6b (and Fig. 10b). The data for uniform heights (both for the aligned and staggered arrays) show a satisfying collapse onto the fully-rough line (black solid line), which is a strong indication that the turbulence statistics are self-similar across cases and that the boundary layer developing above the urban canopies follows the classical scaling laws.

Fig. 6
figure 6

a Mean streamwise velocity profiles in inner scales for the UH cases. The dotted black line represents a smooth wall. b Diagnostic plot for the same cases. The black solid line represents the fully-rough regime (Castro et al. 2013), while the grey line is the smooth-wall limit (Alfredsson et al. 2011)

3.1.3 The Roughness Sublayer

The top of the RSL is considered to be the point where all the effects of the individual roughness elements on the flow cease (Cheng et al. 2007; Kent et al. 2017). It follows that the flow inside the RSL is not spatially homogeneous. The nine vertical profiles taken over the surface are presented in Fig. 7a, b for the UHA and UHS cases, respectively. These are shown to converge unexpectedly at about 1.15\(h_{avg}\) for both cases. The uniform-height surface results for this tall-element canopy differ from the small-cube cases analyzed by Cheng and Castro (2002), where the RSL was found to be much deeper (\(\approx 2h_{avg}\)). Despite the comparable Reynolds number, the RSL depth is much shallower here, closer to that found by Florens et al. (2013). Given the similar conditions between the current data and those of Cheng and Castro (2002), we attribute this difference to the discrepancy in the packing density and the much deeper canopy layer. The current results also vary significantly from the commonly cited range of 2–5\(h_{avg}\) (Kent et al. 2017; Flack et al. 2007; Roth 2000) due to the combination of high-aspect-ratio elements and a dense canopy. Furthermore, Fig. 7a, b shows a clear difference between cases, with Fig. 7b showing similarity and collapse between all profiles. The flow inside the canopy behaves as predicted by the literature, similar to that of a densely forested canopy, with a region of severe velocity deficit up to the elements’ average height; this is accordance with Fig. 9.7 in Sect. 9.7.3 from Stull (1988). These results are qualitatively in agreement with those in Nepf (2012) for vegetation canopies for which an inflection point in the velocity profiles appears just above the roughness height for a dense canopy (i.e. \(\lambda _f>0.23\)). The tightly packed roughness generates a strong shear layer at the top of the elements. Cheng et al. (2007) argues that this region is the main contributor to the magnitude of \(z_{0}\), as demonstrated in the 20-mm cube arrays. Figure 7a, however, highlights the effect of aligned street canyons. Profiles 1, 4, and 7 exhibit flow penetration into the canopy where the street canyons are aligned. The rest of the profiles follow a trend very similar to that in Fig. 7b. This demonstrates that, for \(\lambda _{p} = 0.44\), the flow cannot penetrate deeply into the canopy, unless the streets are aligned. Even so, velocity profiles 1, 4, and 7 only begin to have a significant increase at around 0.5\(h_{avg}\), indicating that at this value of \(\lambda _{p}\), even in aligned street canyons, the penetration into the canopy is limited beyond the depth of 40 mm. A further study could be conducted on varied-canopy-depth models to verify how deep the flow can penetrate in aligned surfaces with varying h at different packing densities. In summary, uniform tall arrays, concerning the RSL, behave similarly to forested canopies, showing an inflection point in correspondence with the top of the canopy. Lastly, for large values of \(\lambda _{p}\) and \(h/\delta \), even when the elements are aligned, the flow cannot penetrate significantly within the canopy.

Fig. 7
figure 7

b, d RSL depth shown with vertical streamwise velocity profiles for uniform-height models. a, c Positions of profiles shown over a module. a, b UHA; c, d UHS. The legends of b, d refer to the profile numbers in a, c

3.1.4 The Inertial Sublayer

In the ISL, the wall-normal variation of shear stress may be neglected, hence the constant-flux-region denomination (Kent et al. 2017). Here, we defined the ISL as the region where the vertical variation of the spatially-averaged profiles of Reynolds shear stresses is below ± 10%, as done for quasi-constant shear stress layer in Cheng et al. (2007). The profiles are reported in Fig. 8a, c. The base of the ISL is assumed to be the top of the RSL found in Sect. 3.1.3. In contrast to Cheng and Castro (2002) and Cheng et al. (2007), where the determination of an ISL was found to be challenging, here an ISL of significant depth can be observed in both cases. Jiménez (2004) reported \(\delta /h\approx 80\) as the lower limit for a canonical logarithmic region to exist. For the uniform-height model, \(\delta /h\approx 3\), well below the limit where an ISL should form. Arguably, in Fig. 8 an ISL forms for these cases. Cheng et al. (2007) also argued that the ISL depth is not constant, which is also demonstrated here. As discussed in Sect. 3.1.1, it is likely that a deep ISL forms due to the large value of \(\lambda _{p}\). It is intuitive to imagine that, in the limit of \(\lambda _{p}\rightarrow \infty \), a new raised smooth wall surface (at \(z=h\)) forms, and recovers the canonical turbulent boundary-layer structure. Additionally, it was discussed in Sect. 3.1.2 how both the mean streamwise velocity profiles and the turbulent statistics are found to conform to the canonical scaling for turbulent boundary layers, and how the roughness parameters are invariant with the streamwise location once the flow is fully-developed. These findings strengthen the argument that, despite the relative roughness height considered, an ISL does indeed develop over these urban surfaces.

Fig. 8
figure 8

a, c The CSL shown with a shear-stress scatter plot and spatially-averaged shear-stress profile (black line). Different marker colours represent different profile locations within a roughness module. Black dotted lines show the CSL boundaries determined by the ± 10% from the averaged CSL value. The blue line indicates \(\overline{u'w'}\) determined by the average CSL value, whilst the red line uses best-fit extrapolation to the value of \(h_{avg}\) to determine \(\overline{u'w'}\) [as in Florens et al. (2013)]. b, d Equation 1 is rearranged to the form of \(z=z_{0}e^{k{U}/{u}_{*}}+d \), where \(u_*\) is calculated using \(\overline{u'w'}\). Using this rearranged form, \(z_0\) represents the gradient and d is the axis intercept. A least-mean-squares fit is used to extrapolate onto the axis within the CSL boundaries. Case UHA for a and b, and case UHS for c and d

3.1.5 The Aerodynamic Parameters

The method used to calculate d and \(z_{0}\) uses a logarithmic-law fitting. A common problem faced with the uncertainty in the fitting of the logarithmic law is that there are three free parameters: \(u_{*}\), d and, \(z_{0}\) (Castro 2007; Segalini et al. 2013). More generally, while other unknowns are given by the wake parameter, \(\varPi \), and by not considering the von Kármán ‘constant’ a constant (Castro et al. 2017), this is not relevant for the current discussion. To reduce the uncertainty in the fitting procedure, a common method involves first fixing the value of \(u_{*}\),and is performed here with two methods: (i) by using the Reynolds shear -stress value in the ISL to determine \(u_{*}\), and (ii) by measuring the drag directly by use of pressure-tapped elements. In the first method, the horizontally-averaged vertical streamwise velocity profile is used in either of two ways. First, one can simply compute the average value from all the points within the ISL; second, a best-fit extrapolation onto the height \(z = d\) or \(z = h\) can be used (Cheng and Castro 2002; Cheng et al. 2007; Florens et al. 2013). An independent method to compute the value of \(u_{*}\) is with the use of a pressure-instrumented element, as described in Cheng and Castro (2002), with the results summarized in Table 3. Once the friction velocity has been determined, the aerodynamic parameters are computed by best fitting the velocity profiles considering the von Kármán constant \(\kappa = 0.4\), which is within the limits suggested by Marusic et al. (2013) for high Reynolds number flows. The discrepancy between the friction velocities calculated via averaged and extrapolated ISL is below 4% (Fig. 8b, d). A further comparison extrapolating to \(z = h\) (Florens et al. 2013) gives results more closely resembling the pressure-tapped results (i.e. \(\approx \) 1.5% and \(\approx \) 15% for d and \(z_{0}\) respectively), suggesting that the extrapolation technique does yield to more robust results.

3.2 Effect of Varying Element Height

3.2.1 Depth of the Boundary Layer

In the varied-height models, similar to the previous case, an initial investigation was carried out to determine the fetch at which the boundary layer is fully developed and this is found to be at \(x = 4470~\hbox {mm}\) (\(\approx 9\delta \)), as in Fig. 9. Although the value of the \(h_{avg}\) is the same in both models, the boundary-layer depth is approximately double in the varied-height model (\(\delta _{UHA}/h_{avg} = 3\) versus \(\delta _{VHA}/h_{avg} = 6.25\)). It is likely that the increased value of \(\sigma _{h}\) also increases the drag generated by the surface, hence creating deeper boundary layers.

Fig. 9
figure 9

Boundary-layer depth found from vertical streamwise velocity profiles taken at different fetch over the models a VHA and b VHS; x-units are in mm

3.2.2 Mean and Fluctuating Streamwise Velocity Profiles

As for the case of uniform-height previously discussed, mean and fluctuating velocity profiles are considered next. Mean streamwise velocity profiles for flow over varied-height cases are shown in Fig. 10a. As for the previous case, a region of near-linear behaviour is noticeable, however, both the degree of collapse onto the smooth-wall dashed line and the extent of the linear regime are reduced from the previous cases in Fig. 6a. While both these aspects are discussed in the following sections, these seem to indicate that the methodology for the evaluation of the roughness parameters (see Sect. 3.2.4) has a higher uncertainty for surfaces with height heterogeneity and that the roughness sublayer has grown at the expense of the ISL (see Sect. 3.2.3), respectively. It is still important to highlight that, even for this cases, data from different streamwise stations yield converged roughness parameters, which is a strong evidence of the fully-developed nature of the flow at the measuring stations.

Fig. 10
figure 10

a Mean streamwise velocity profiles in inner scales for the VH cases. The dotted black line represents a smooth wall. b Diagnostic plot for the same cases. The black solid line represents the fully-rough regime (Castro et al. 2013), while the grey line is the smooth-wall limit (Alfredsson et al. 2011)

The diagnostic plot for the streamwise fluctuations is shown in Fig. 10b. Both cases for varied heights are found to be self-similar, however interestingly, these sit above the fully-rough trend line. This can be interpreted as an indication of the enhanced turbulent mixing due to the heterogeneous surface morphology, an aspect further discussed in Sect. 3.2.3. This discrepancy with previous data on fully-rough walls reported in Castro et al. (2013) can also be attributed to the possibility that this fully-rough asymptote within the diagnostic plot is effected by several other roughness parameters, as highlighted originally in Castro et al. (2013). Two parameters of particular relevance here are the much higher relative roughness height, \(\delta /h\), and the standard deviation in elements height, \(\sigma _h\), which distinguish this work from data in Castro et al. (2013).

3.2.3 The Roughness Sublayer

In the varied-height experiment, 81 profiles were taken per repeating \(3\times 3\) modules (Fig. 2), and are shown in Fig. 11a, b for the two configurations. There is a clear collapse of mean vertical profiles, indicating the existence of a limited RSL depth (Fig. 11a, b). The RSL depth is the same in both VHA and VHS configurations, which is of interest; it seems that the dominant feature in determining the depth of the RSL is the tallest element in the \(3\times 3\) module, as a collapse is found for \(z/h_{avg}\approx 2.5\) (\(z/h_{max}\approx 1\)). The evident collapse of vertical profiles in the RSL contradicts Cheng and Castro (2002), Cheng et al. (2007), and Jiménez (2004), who predicted an RSL region that expands and ‘squeezes’ the ISL for rough walls with significantly tall and heterogeneous elements. For the VHA surface (Fig. 11a), not all profiles collapse at the same height. Several vertical profiles converge near the average height of the elements (i.e., \(z/h_{avg} = 1\)). The remaining profiles converge just above the maximum-height element (\(z/h_{avg} = 2.5\) or \(z/h_{max}\approx 1\)). This trend possibly occurs for the same reasons discussed for the uniform-height cases. Where street canyons are aligned in the wind direction, the flow penetrates deeper into the canopy, allowing for higher velocities and enhanced mixing. Likewise, in the VHS configuration (Fig. 11b), the flow does not seem to fully converge until above the tallest element (\(z/h_{avg} = 2.85\)). In the varied-height situation, skimming no longer occurs since there is a large spread of velocities below the heights \(h_{avg}\) and \(h_{max}\). As seen in both the VHA and VHS cases, the large spread of velocities within the canopy indicates significant flow penetration likely due to the increased \(\sigma _{h}\). Considering this possibility, cities designed with large \(\sigma _{h}\) could introduce much higher rates of mixing, with positive outcomes on urban air quality and natural ventilation. As argued by Zaki et al. (2011), \(\lambda _{p}\) is no longer a good indicator of the flow regime when a large \(\sigma _{h}\) is introduced, which is counter to the studies conducted by Sharma and García-Mayoral (2020). The argument made here suggests that the surface parameters used to describe flow over uniform-height models are no longer sufficient to characterize canopies with significant height variations.

Fig. 11
figure 11

RSL shown with vertical streamwise velocity profiles over the varied-height models. a VHA and b VHS. Horizontally-averaged vertical streamwise velocity profiles of different experimental runs. c VHA d VHS; x-units are in mm

Figure 11c, d shows horizontally-averaged profiles within the RSL for the different experimental runs in the varied-height model. The spatially-averaged profiles taken in different areas of the array (across the wind-tunnel span) collapse reasonably well onto each other, but there are some discrepancies around the region \(1.5<z/h<2.5\). The profiles measured across the y-planes collapse better than the profiles taken over the \(3\times 3\) module. Possibly, it is more efficient and easier to obtain a representative averaged profile from taking measurements across the width of the model rather than over a single repeating module. Another reason for this better collapse could be that the average height of the elements in front of and behind the y-planes does not always match \(h_{avg}\)—see further discussion on this topic in Cheng and Castro (2002). However, the clear collapse in y-plane profiles indicates that \(h_{avg}\) may not be the dominant length scale in models with large \(\sigma _{h}\). Comparing uniform- to varied-height canopies across Figs. 7 and 11, a prominent difference appears. For the varied-height models, the flow velocities vary greatly within the canopy, in contrast, they remain close to zero until the very top of the canopy for the uniform-height experiments. This is likely due to the dynamics of the ‘skimming flow’ regime (Grimmond and Oke 1999). Even below \(h_{avg}\), there are much higher velocities in the varied-height canopy compared with the uniform-height. The spread of velocity profiles in the varied-height model suggests that reasonable mixing can occur in the near-wall region for a tall and dense canopy (large \(\lambda _{p}\) and \(h_{avg}\)), providing that the magnitude of \(\sigma _{h}\) is significant.

3.2.4 The Inertial Sublayer

Several thresholds have been used in the literature for defining an ISL, where ± 5, ± 10, and ± 20% variation are all examined by Cheng et al. (2007). Kanda et al. (2013) used a region for logarithmic-law fitting, which is a function of the tallest and average roughness elements’ height. Sharma and García-Mayoral (2020) noted how the ISL depth in their work is a function of the spacing between the roughness elements. These ISL fitting methods, however controversial, can prove very useful in accurately estimating aerodynamic parameters. However, a convincing relationship between the roughness geometry and the upper limits of the ISL is still elusive, particularly for tall canopies. The results discussed in Sect. 3.2.2 are somewhat inconclusive on the existence of a self-similar logarithmic region developing over tall rough surfaces of varying height as the roughness sublayer gradually replaces the ISL/CSL. This is visible both in the logarithmic portion of the wind profiles, which is not as extensive, or convincing as before, and by the higher-order quantities that are found to deviate from the experimental trend within the diagnostic plot (Fig. 10). This is in stark contrast with the previous cases of uniform height, which are shown to behave according to the theory in Fig. 6. However, our analysis also suggested that the roughness parameters (calculated from the logarithmic region) show independence from both the boundary-layer height and the spanwise location of the measurements; this would offer support to the hypothesis of a true ISL developing over these surfaces despite their low \(\delta /h\) ratio. For these reasons, we prefer referring to a CSL rather than an ISL for the cases of varying height in the remainder of this section. We believe that the knowledge of the existence of a CSL would still be beneficial in the modelling of tall heterogeneous canopies.

Fig. 12
figure 12

a, c The CSL shown with a shear-stress scatter plot and spatially-averaged shear-stress profile (black line). Different marker colours represent different profile locations within a roughness module. Black dotted lines show the CSL boundaries determined by the ± 10% from the averaged CSL value. The blue line indicates \(\overline{u'w'}\) determined by the average CSL value, whilst the red line uses best-fit extrapolation to the value of \(h_{avg}\) to determine \(\overline{u'w'}\) [as in Florens et al. (2013)]. b, d Equation 1 is rearranged to the form of \(z=z_{0}e^{k{U}/{u}_{*}}+d\), where \(u_*\) is calculated using \(\overline{u'w'}\). Using this rearranged form, \(z_0\) represents the gradient and d is the axis intercept. A least-mean-squares fit is used to extrapolate onto the axis within the CSL boundaries. Case VHA for a and b, and case VHS for c and d

The depth of the CSL region is shown in Fig. 12. Previous studies predicted that an ISL region would vanish as the RSL rises through the boundary layer due to increasingly rough surfaces (Cheng et al. 2007; Cheng and Castro 2002; Rotach 1999; Jiménez 2004; Hagishima et al. 2009). The clear collapse of the RSL in Fig. 11 partly contradicts this theory. Furthermore, the ± 10% variation definition does allow for a CSL with significant depth to be singled out, as seen in Fig. 12. This may be an indication that the CSL does, indeed, exist. It is important to point out, however, that the Reynolds shear stresses never reach a full plateau, but show a small—yet visible—gradient across the CSL. This is possibly related to the fact that the boundary layer in the cases with varied height occupies a significant height of the wind tunnel. The acceleration parameter is fairly small (9.48 \(\times \ 10^{-8}\)) indicating a zero-pressure gradient, however, it must be acknowledged that it is nearly double the same quantity in the uniform-height case. A CSL within the ± 10% variation does form even though \(\delta /h= 6.25\) is still well below the limits established by Jiménez (2004) for the ISL. Florens et al. (2013) argued that this ratio should be lowered, and additionally Cheng et al. (2007) demonstrate that the ratio \(\delta /h\) should not be the only criterion for the existence of the ISL. Furthermore, Leonardi and Castro (2010) argued that the Reynolds number previously thought to be necessary for the ISL development may be much lower than expected. The finding from the current work, therefore, offers some support to the idea that Jiménez ’s (2004) ratio is perhaps too stringent, certainly for a CSL. The \(\pm 10\%\) horizontal variation may not, however, be the best range for determining the CSL depth. The top limit of the CSL can be hard to determine, and it is here tempting to soften the criterion to include the region where the gradient of the slope is more noticeable (i.e., \(2.5<z/h<5\) in Fig. 12c). Indeed, in the uniform-height experiment discussed earlier, the ±10% horizontal variation is quite liberal, and a more constrained ±5% could be used. The depth of the CSL could perhaps be better determined by the change in slope in the spatially-averaged Reynolds shear stresses, in a procedure similar to that in Kanda et al. (2013). It is important to add a note to this discussion to clarify that the terms CSL and ISL are, at times, used instinctively in the literature, while they indicate different entities herein, as previously discussed. A summary of the depths of the CSL and RSL for both cases of heterogeneous heights can be seen in Table 3. The heights of the RSL, CSL, and boundary layer seem to be strongly dependent on the height of the tallest element within the surfaces. Even with appropriate normalization (e.g., by \(h_{avg}\) or \(h_{max}\)), there is no apparent universal scaling between results from varied- and uniform-height experiments. The significant change in the boundary layer, RSL, and CSL depths across the two cases seem to be uncorrelated. Thus, VHA and VHS surfaces could only be compared to models of similarly sized standard deviation and average height.

3.2.5 The Aerodynamic Parameters

The aerodynamic parameters and friction scaling were determined as in Sect. 3.1.5, but direct measurement of the friction velocity by an instrumented element were not available for these cases for reasons discussed in Sect. 2.3. Here we follow the same procedure explained in Sect. 3.1.5, i.e., extrapolating to \(z = h\), with the results reported in Fig. 12. The aerodynamic parameters calculated herein are compared to the predictions of the Kanda et al. (2013) and the Macdonald (2000) morphometric methods, which considered elements with varied height, to assess their performance. There is a large discrepancy between the varied- and uniform-height results, where d is found to increase 1.5 times from uniform- to varied-height, while \(z_{0}\) increases nearly 3.5 times. For the uniform canopy, the Kanda et al. (2013) method predicted reasonable d values, while the Macdonald (2000) performed worse (\(\approx \) 5% and 30%, respectively). Neither produced a value for normalized \(z_0\) that is close to that obtained herein (with both methods predicting a roughness length almost an order of magnitude larger). For the varied-height experiment, the extrapolated results from the shear-stress plot produce a value of d above \(h_{max}\), while the average shear-stress method results in values below \(h_{max}\) (see Fig. 12). Both the Kanda et al. (2013) and Millward-Hopkins et al. (2012) methods predict a value of d between h and \(h_{max}\), which differ \(\approx 36\%\) and \(\approx 43\%\), respectively from those extrapolated here from the average shear. This discrepancy may be due to the fact that Tokyo (heavily used in Kanda et al. (2013)) does not have many super-tall structures, such as Hong Kong (largely due to the seismic restrictions in the former). The lack of similarity of results continues when comparing \(z_{0}\). The \(z_{0}\) values calculated with the two methods were roughly twice as large those evaluated here.

Our results strengthen the argument by Kanda et al. (2013) that new parameters are necessary to accurately characterize an urban rough wall, particularly when the elements are tall and heterogeneous in height if one wants to accurately estimate the aerodynamic parameters. Jiang et al. (2008) and Zaki et al. (2011) suggest that the magnitude of d increases with \(\sigma _{h}\), which is in line with results presented here. Zaki et al. (2011) also demonstrated that the values of d and \(z_{0}\) almost double for cases with the same average height but significant \(\sigma _{h}\) variation. Xie et al. (2008) also acknowledge that the tallest elements have a disproportionate contribution to the drag across a surface. Kanda et al. (2013) further suggested that both \(\lambda _{p}\) and \(\sigma _{h}\) values have an effect on d and \(z_{0}\). A lack of similarity between results of uniform- and varying-height strengthens the above arguments, particularly the fact that a large \(\sigma _{h}\) produces significant changes in the flow field. Kanda et al. (2013) considered \(h_{max}\) as the upper limit of d, and also advocate the use of the standard deviation of element height (\(\sigma _{h}\)) as an important measure when calculating the aerodynamic parameters, which is corroborated by the findings here. This seems to support the argument that the sole use of \(\sigma _{h}\) and \(h_{avg}\) when designing an urban model is not representative of the flow dynamics in a real city as suggested by Kanda et al. (2013). Since there are usually a few very tall elements amongst many shorter elements, other important parameters to consider would be the distribution and ratios of tall to short elements.

4 Conclusions

Wind-tunnel experiments were conducted at the University of Surrey on four dense (\(\lambda _p=0.44\)) and tall (\(\delta /h_{avg}\approx 3\)) urban arrays; two canopies were of uniform height and two of varied height, while the average height was kept fixed at \(h_{avg}\) = 80 mm in both cases. All canopies aimed to represent idealized modern cities. The experiments examined the differences in the flow features across canopies with homogeneous and heterogeneous heights by means of LDA and direct-drag measurements. All cases examined were fully-developed, fully-rough surfaces in zero-pressure gradients.

In the uniform-height experiments, the surfaces with a large value of \(\lambda _{p}\) inhibits deep penetration of the flow into the canopy, thus hindering mixing at street level. For \(\lambda _{p} = 0.44\), a skimming-flow regime seemingly occurs, and the rough surface starts to recover smooth-like properties above the elements, producing an inertial sublayer region; this is in contrast to Cheng and Castro (2002) and Cheng et al. (2007) who questioned the existence of this region as the roughness influence grows. The RSL depth is found to extend approximately up to \(z = 1.15h\), which is much shallower than the typically expected 2–5h (Flack et al. 2007). Conventional turbulent scalings laws are found to be applicable to both mean streamwise velocity profiles and turbulent quantities despite the severity of the roughness.

For heterogeneous-height canopies, the usefulness of \(\lambda _{p}\) in describing the wall properties becomes more questionable. Despite the same \(h_{avg}\), the boundary layer grows to almost double the thickness of the cases with uniform height, highlighting the significance of the parameters \(\sigma _{h}\) and \(h_{max}\) in dictating the flow features. Surprisingly, even with such a heterogeneous canopy, a clear collapse of vertical profiles of Reynolds shear stress is observed, forming a coherent RSL extending to just above the height of the tallest element (\(z/h_{avg}=2.85\) and \(z/h_{max}\approx 1\)), which is much shallower than anticipated. Moreover, more significant flow penetration is observed within the canopy when compared with the uniform-height array. This is responsible for an enhanced turbulent mixing resulting in velocity fluctuations that are higher than previously reported in surfaces with homogeneous heights throughout the boundary-layer depth (Castro et al. 2013). These findings strengthen the need to include information regarding \(\sigma _{h}\) and \(h_{max}\) when describing the flow over tall and heterogeneous canopies, supporting the conclusions highlighted in Kanda et al. (2013). A comparison of the aerodynamic parameters for the cases with heterogeneous height considered herein with existing morphometric methods has highlighted the inaccuracies of these methods in the case of tall canopies with significant variation in height between elements. Further work is therefore needed on this topic.