1 Introduction
Canopy flows are ubiquitous in both natural and artificial settings. Although they have mostly been studied in the framework of flows through crops and forests (Finnigan Reference Finnigan2000; Belcher, Harman & Finnigan Reference Belcher, Harman and Finnigan2012; Nepf Reference Nepf2012), they are also relevant to flows over engineered surfaces, such as pin fins for heat transfer or piezoelectric filaments for energy harvesting (Bejan & Morega Reference Bejan and Morega1993; McCloskey, Mosher & Henderson Reference McCloskey, Mosher and Henderson2017). While the latter usually have small to moderate Reynolds numbers, vegetation canopies typically have much larger ones. The study of turbulent flows over canopies has wide-ranging applications, including reducing crop loss (de Langre Reference de Langre2008), energy harvesting (McGarry & Knight Reference McGarry and Knight2011; McCloskey et al. Reference McCloskey, Mosher and Henderson2017; Elahi, Eugeni & Gaudenzi Reference Elahi, Eugeni and Gaudenzi2018) and improving heat transfer (Fazu & Schwerdtfeger Reference Fazu and Schwerdtfeger1989; Bejan & Morega Reference Bejan and Morega1993).
On the basis of the geometry and spacing of the canopy elements, a canopy can be classified as dense, sparse or transitional (Nepf Reference Nepf2012). In the dense limit, the canopy elements are in close proximity to each other and turbulence is essentially not able to penetrate within the canopy layer. In the sparse limit, the spacing between canopy elements is large and the turbulent eddies can penetrate the full depth of the canopy. An intermediate or transitional regime lies between these two limits. Turbulent flows in the dense regime, reviewed by Finnigan (Reference Finnigan2000) and Nepf (Reference Nepf2012), are characterised by the formation of Kelvin–Helmholtz-like, or mixing-layer, instabilities, originating from the inflection point at the canopy tips (Raupach, Finnigan & Brunet Reference Raupach, Finnigan and Brunet1996). As the sparsity of the canopies is increased, the importance of the Kelvin–Helmholtz-like instability decreases (Poggi et al. Reference Poggi, Porporato, Ridolfi, Albertson and Katul2004; Huang, Cassiani & Albertson Reference Huang, Cassiani and Albertson2009; Pietri et al. Reference Pietri, Petroff, Amielh and Anselmet2009). Eventually, the flow would resemble that over a smooth wall, albeit perturbed by the discrete presence of the individual canopy elements (Finnigan Reference Finnigan2000). The separation between these regimes is still somewhat unclear. Nepf (Reference Nepf2012) proposed an approximate classification of the canopy regime based on the roughness frontal density, 
$\unicode[STIX]{x1D706}_{f}$. Nepf observed that canopies are dense for 
$\unicode[STIX]{x1D706}_{f}\gg 0.1$, sparse for 
$\unicode[STIX]{x1D706}_{f}\ll 0.1$, and intermediate for 
$\unicode[STIX]{x1D706}_{f}\approx 0.1$. However, in addition to the geometric parameter 
$\unicode[STIX]{x1D706}_{f}$, the length scales of the flow should also be considered when determining the regime of the canopy. The length scales in a turbulent flow may be much larger than the element spacing at a particular Reynolds number, so that the turbulent eddies are precluded from penetrating within the canopy. As the Reynolds number is increased, however, the turbulent length scales will eventually become comparable to the element spacing and allow the turbulent eddies to penetrate within the canopy efficiently.
In the present work, we study flows within and above sparse canopies using direct numerical simulation (DNS). While this allows for the full resolution of turbulence, it restricts our simulations to moderate friction Reynolds numbers, 
$Re_{\unicode[STIX]{x1D70F}}\approx 520{-}1000$, and element heights, 
$h^{+}\approx 110{-}200$. While these heights would be directly applicable to some of the engineered canopies mentioned above, they are much smaller than those typical of vegetation canopies, 
$h^{+}\approx 10^{4}{-}10^{6}$ (e.g. Green, Grace & Hutchings Reference Green, Grace and Hutchings1995; Novak et al. Reference Novak, Warland, Orchansky, Ketler and Green2000; Zhu et al. Reference Zhu, Van Hout, Luznik, Kang, Katz and Meneveau2006), although comparable to some laboratory experiments, 
$h^{+}\approx 400{-}800$ (Raupach, Hughes & Cleugh Reference Raupach, Hughes and Cleugh2006; Böhm et al. Reference Böhm, Finnigan, Raupach and Hughes2013). In any event, we provide evidence in § 2.3 of the scaling of the canopy-flow dynamics with the Reynolds number, which would make our conclusions of relevance for canopies at larger Reynolds numbers as well. The present canopies have low roughness densities 
$\unicode[STIX]{x1D706}_{f}\lesssim 0.1$, with element spacings large enough to limit their interaction with the near-wall turbulence dynamics. Owing to the sparse nature of these canopies, we would expect the flow within them to be dominated by the footprint of the canopy elements, rather than by a mixing-layer instability.
Conventionally, a homogeneous drag is used to represent the effect of canopies (Dupont & Brunet Reference Dupont and Brunet2008; Finnigan, Shaw & Patton Reference Finnigan, Shaw and Patton2009; Huang et al. Reference Huang, Cassiani and Albertson2009; Bailey & Stoll Reference Bailey and Stoll2016). This approach would only be strictly valid to represent very closely packed canopies, where the element spacing is much smaller than any length scale in the overlying flow, and even small flow structures perceive the canopy elements as acting collectively (Zampogna & Bottaro Reference Zampogna and Bottaro2016). Using a homogeneous drag to capture the effects of sparser canopies tends to overdamp turbulent fluctuations within the canopies (Yue et al. Reference Yue, Parlange, Meneveau, Zhu, van Hout and Katz2007; Bailey & Stoll Reference Bailey and Stoll2013). This is typically attributed to the inability of homogenised models to capture the element-induced flow, and the lack of representation of the gaps between the canopy elements, where the fluctuations would not experience any damping (Bailey & Stoll Reference Bailey and Stoll2013).
In the present work we separate the effect of the element-induced coherent flow from the incoherent background turbulence, and focus mainly on the properties of the latter. We study different element spacings and geometries. We propose a scaling that suggests that the dynamics of the background turbulence within sparse canopies is mainly governed by their effect on the mean velocity, rather than by the direct interaction of the canopy elements with the flow. Based on this scaling, we propose that the effect on the background turbulence is represented better by a drag acting on the mean flow alone than by a homogeneous drag. Partial results from some of the simulations can be found in Sharma & García-Mayoral (Reference Sharma and García-Mayoral2018a,Reference Sharma and García-Mayoralb).
The paper is organised as follows. The numerical methods used and the canopy geometries simulated are described in § 2. The results of the canopy-resolving simulations and the scaling of turbulent fluctuations are discussed in § 3. The results obtained from simulations that substitute the canopy by a drag force are discussed in § 4. The conclusions are presented in § 5.
Figure 1. Schematic of the numerical domain. An instantaneous field of the fluctuating streamwise velocity, scaled, from case TP1 is shown in three orthogonal planes.
2 Numerical simulations
We conduct DNS of an open channel with canopy elements protruding from the wall. The streamwise, wall-normal and spanwise coordinates are 
$x$, 
$y$ and 
$z$, respectively, and the associated velocities are 
$u$, 
$v$ and 
$w$. The size of the simulation domain is 
$2\unicode[STIX]{x03C0}\unicode[STIX]{x1D6FF}\times \unicode[STIX]{x1D6FF}\times \unicode[STIX]{x03C0}\unicode[STIX]{x1D6FF}$, with the channel height 
$\unicode[STIX]{x1D6FF}=1$. This box size has been shown to be adequate to capture one-point statistics up to the channel height for the friction Reynolds numbers used in the present study (Lozano-Durán & Jiménez Reference Lozano-Durán and Jiménez2014). A schematic representation of the numerical domain is shown in figure 1. The domain is periodic in the 
$x$ and 
$z$ directions. No-slip and impermeability conditions are applied at the bottom boundary, 
$y=0$, and free slip and impermeability at the top, 
$y=\unicode[STIX]{x1D6FF}$. It is shown in § 3 that the height of the roughness sublayer for the canopies studied here extends to only half of the domain height, so the top boundary of the channel does not interfere with the canopy flow. The flow is incompressible, with the density set to unity. All simulations are run at a constant mass flow rate, with the viscosity adjusted to obtain the desired friction Reynolds number based on the total stress. Most simulations are conducted at a friction Reynolds number 
$Re_{\unicode[STIX]{x1D70F}}=u_{\unicode[STIX]{x1D70F}}\unicode[STIX]{x1D6FF}/\unicode[STIX]{x1D708}\approx 520$, and a few at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$.
The numerical method used to solve the three-dimensional Navier–Stokes equations is adapted from Fairhall & García-Mayoral (Reference Fairhall and García-Mayoral2018). A Fourier spectral discretisation is used in the streamwise and spanwise directions. The wall-normal direction is discretised using a second-order centred difference scheme on a staggered grid. For the simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$, the grid in the wall-normal direction is stretched to give a resolution 
$\unicode[STIX]{x0394}y_{min}^{+}\approx 0.2$ at the wall, stretching to 
$\unicode[STIX]{x0394}y_{max}^{+}\approx 2$ at the top of the domain. For the simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$, wall-normal resolutions of 
$\unicode[STIX]{x0394}y_{min}^{+}\approx 0.35$ and 
$\unicode[STIX]{x0394}y_{max}^{+}\approx 5.5$ are used, as dissipation occurs at larger scales near the centre of the channel at larger Reynolds numbers (Jiménez Reference Jiménez2012). The wall-parallel resolutions for the different cases are given in table 1. The time advancement is carried out using a three-step Runge–Kutta method with a fractional step, pressure correction method that enforces continuity (Le & Moin Reference Le and Moin1991):
$$\begin{eqnarray}\displaystyle \displaystyle \left[\unicode[STIX]{x1D644}-\unicode[STIX]{x0394}t\frac{\unicode[STIX]{x1D6FD}_{k}}{Re}L\right]\boldsymbol{u}_{k}^{n} & = & \displaystyle \boldsymbol{u}_{k-1}^{n}+\unicode[STIX]{x0394}t\left[\frac{\unicode[STIX]{x1D6FC}_{k}}{Re}L\boldsymbol{u}_{k-1}^{n}-\unicode[STIX]{x1D6FE}_{k}N(\boldsymbol{u}_{k-1}^{n})\right.\nonumber\\ \displaystyle & & \displaystyle \left.-\,\unicode[STIX]{x1D701}_{k}N(\boldsymbol{u}_{k-2}^{n})-(\unicode[STIX]{x1D6FC}_{k}+\unicode[STIX]{x1D6FD}_{k})G(p_{k}^{n})\vphantom{\frac{\unicode[STIX]{x1D6FC}_{k}}{Re}}\right]\!,\quad k\in [1,3],\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle DG(\unicode[STIX]{x1D719}_{k}^{n})=\frac{1}{(\unicode[STIX]{x1D6FC}_{k}+\unicode[STIX]{x1D6FD}_{k})\unicode[STIX]{x0394}t}D(\boldsymbol{u}_{k}^{n}), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{u}_{k+1}^{n}=\boldsymbol{u}_{k}^{n}-(\unicode[STIX]{x1D6FC}_{k}+\unicode[STIX]{x1D6FD}_{k})\unicode[STIX]{x0394}t\,G(\unicode[STIX]{x1D719}_{k}^{n}), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle p_{k+1}^{n}=p_{k}^{n}+\unicode[STIX]{x1D719}_{k}^{n}. & \displaystyle\end{eqnarray}$$ Here 
$\unicode[STIX]{x1D644}$ is the identity matrix; 
$L$, 
$D$ and 
$G$ are the Laplacian, divergence and gradient operators, respectively; 
$N$ is the dealiased advective term; 
$\unicode[STIX]{x1D6FC}_{k}$, 
$\unicode[STIX]{x1D6FD}_{k}$, 
$\unicode[STIX]{x1D6FE}_{k}$ and 
$\unicode[STIX]{x1D701}_{k}$ are the Runge–Kutta coefficients for substep 
$k$ from Le & Moin (Reference Le and Moin1991); and 
$\unicode[STIX]{x0394}t$ is the time step.
Table 1. Simulation parameters: 
$N_{x}$ and 
$N_{z}$ are the number of canopy elements in the streamwise and spanwise directions, respectively; 
$u_{\unicode[STIX]{x1D70F}}$ is the friction velocity based on the net drag and scaled with the channel bulk velocity; 
$Re_{\unicode[STIX]{x1D70F}}$ is the friction Reynolds number based on 
$u_{\unicode[STIX]{x1D70F}}$ and 
$\unicode[STIX]{x1D6FF}$; 
$\unicode[STIX]{x1D706}_{f}$ is the roughness frontal density; and 
$\int D^{+}$ is the net canopy drag force scaled with 
$u_{\unicode[STIX]{x1D70F}}$, that is, the proportion of the total drag on the fluid exerted by the canopy elements, with the remainder being the friction at the bottom wall. The grid resolutions in the streamwise and spanwise directions are 
$\unicode[STIX]{x0394}x^{+}$ and 
$\unicode[STIX]{x0394}z^{+}$, respectively.
2.1 Canopy-resolving simulations
We have considered flows over canopies with both permeable and impermeable elements. The geometry of the impermeable canopy elements is resolved using an immersed-boundary method adapted from García-Mayoral & Jiménez (Reference García-Mayoral and Jiménez2011). The simulation parameters for the different cases studied here are summarised in table 1. Case S is an open channel with a smooth-wall floor. The canopy-resolving simulations include two canopy geometries, as portrayed in figure 2, with varying element spacings. The first geometry, denoted by the letter ‘P’, consists of collocated prismatic posts with a square top-view cross-section with sides 
$\ell _{x}^{+}=\ell _{z}^{+}\approx 20$, and height 
$\ell _{s}^{+}\approx 110$. The spacing between the canopy elements in the wall-parallel directions for cases PD, P0, P1 and P2 are 
$L_{x}^{+}=L_{z}^{+}\approx 50$, 100, 200 and 400, respectively. The canopy of case PD has a frontal area density of 
$\unicode[STIX]{x1D706}_{f}\approx 0.88$, which would place it in the dense regime (Nepf Reference Nepf2012). We use this simulation to contrast sparse and dense canopy dynamics. The second geometry, denoted by the letter T, consists of frontally extruded T-shaped canopies, as portrayed in figure 2(b,c), in a collocated arrangement. We consider two element spacings for the T-shaped canopies, with cases T1 and T2 having 
$L_{x}^{+}=L_{z}^{+}\approx 200$ and 400, respectively. The head of these canopy elements has dimensions 
$\ell _{x}^{+}=\ell _{z}^{+}\approx 40$ in the wall-parallel directions. The base of the canopy elements has 
$\ell _{x}^{+}\approx 40$ and 
$\ell _{z}^{+}\approx 20$. The base and the head are 
$\ell _{s}^{+}\approx 80$ and 
$\ell _{h}^{+}\approx 30$ tall, respectively.
We also study how canopies with permeable and impermeable canopy elements affect the surrounding flow. The permeable canopy elements of case TP1 have the same geometry and layout as T1, but the elements are represented by a drag force only applied within them, rather than by immersed boundaries. This method allows some flow to permeate into the canopy elements, as shown in figure 2(c,f), and has been observed to be a suitable model for certain natural canopies (Yue et al. Reference Yue, Parlange, Meneveau, Zhu, van Hout and Katz2007; Yan et al. Reference Yan, Huang, Miao, Cui and Zhang2017). The drag force in case TP1, applied only at the grid points that are within the canopy elements, is of the form 
$C_{dc}u_{i}|u_{i}|$, similar to Yue et al. (Reference Yue, Parlange, Meneveau, Zhu, van Hout and Katz2007), Bailey & Stoll (Reference Bailey and Stoll2013) and Yan et al. (Reference Yan, Huang, Miao, Cui and Zhang2017), where 
$C_{dc}$ is a drag coefficient and 
$u_{i}$ is the instantaneous local velocity in every 
$i$ direction. The value of 
$C_{dc}$ is set such that further increasing its magnitude does not significantly increase the net drag force on the mean flow. This forcing provides a local body force opposing the flow inside the canopy elements, and thus results in a small velocity within the canopy elements. The net mean drag force for this canopy is similar to that of the impermeable canopy, T1, as noted in table 1, in spite of the different character of the canopy elements.
Figure 2. Contours of instantaneous streamwise velocity in planes passing through the centre of a canopy element. Panels (a–c) represent cuts in the 
$z{-}y$ plane, and (d–f) represent cuts in the 
$x{-}y$ plane. Panels are for cuboidal canopy elements from case P1 (a,d); and for T-shaped canopy elements from case T1 (b,e) and case TP1 (c,f). The white lines mark the positions of the canopy elements. The contours are scaled using the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, of each case.
In order to ascertain the effect of the Reynolds number on the results, two additional simulations, 
$\text{P2I}_{Re}$ and 
$\text{P2O}_{Re}$, are conducted at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$. The canopy of 
$\text{P2I}_{Re}$ matches the parameters of the canopy of P2 in inner units, that is, element widths 
$\ell _{x}^{+}=\ell _{z}^{+}\approx 20$, height 
$\ell _{s}^{+}\approx 110$ and element spacings 
$L_{x}^{+}=L_{z}^{+}\approx 400$. The ratio of channel to canopy height for 
$\text{P2I}_{Re}$ is 
$\unicode[STIX]{x1D6FF}/\ell _{s}\approx 10$, whereas for case P2 the ratio is 
$\unicode[STIX]{x1D6FF}/\ell _{s}\approx 5$. The simulation 
$\text{P2I}_{Re}$ is conducted to verify that the channel height used is large enough not to constrain the canopy-layer dynamics. The canopy of 
$\text{P2O}_{Re}$ matches the parameters of the canopy of P2 in outer units, that is, 
$\ell _{x}/\unicode[STIX]{x1D6FF}=\ell _{z}/\unicode[STIX]{x1D6FF}\approx 0.04$, height 
$\ell _{s}/\unicode[STIX]{x1D6FF}\approx 0.2$ and element spacings 
$L_{x}/\unicode[STIX]{x1D6FF}=L_{z}/\unicode[STIX]{x1D6FF}\approx 0.8$. This simulation is conducted to assess Reynolds-number effects for a fixed canopy geometry.
The roughness densities of the canopies are given in table 1. All the canopies studied lie within the sparse to transitional regime empirically demarcated by Nepf (Reference Nepf2012), except that of case PD, which lies in the dense regime. The spanwise spacings between the sparse canopy elements are 
$L_{z}^{+}\gtrsim 100$, which is comparable to, or larger than, the width of near-wall streaks, 
$\unicode[STIX]{x1D706}_{z}^{+}\approx 100$ (Kline et al. Reference Kline, Reynolds, Schraub and Runstadler1967). This implies that the canopies should be sparse from the point of view of the near-wall turbulent fluctuations as well.
Figure 3. Drag coefficients, 
$C_{dh}=D/U^{2}$, obtained from T-shaped canopies: red ——, case T1; violet ——, case TP1; blue ——, case T2; red – – –, cases T1-H/H0; and blue – – –, case T2-H0. The inset provides a magnified view of the drag coefficients for cases T2 and T2-H0.
2.2 Drag force representations
In order to explore the canopy-flow dynamics, we also conduct simulations where the canopy is replaced by some drag force, which does not resolve the geometry of the canopy elements. Sparse canopies consisting of bluff elements, such as those in the present work, are generally better characterised by a quadratic, form drag (Coceal, Thomas & Belcher Reference Coceal, Thomas and Belcher2008), whereas for canopies with slender, filamentous elements, where viscous effects can dominate, a drag force proportional to the velocity would be more appropriate (Tanino & Nepf Reference Tanino and Nepf2008; Sharma & García-Mayoral Reference Sharma and García-Mayoral2019). For complex natural canopies with foliage, which can have a range of element scales, both form and viscous drags can be important (Finnigan Reference Finnigan2000). For the canopies studied here, we find that the drag is essentially quadratic, and thus replace the canopy elements by a quadratic drag force. Note, however, that whether the drag was linear, quadratic, or otherwise is of no consequence to the conclusions that we derive.
The drag coefficient, 
$C_{dh}$, is obtained by approximating the canopy drag force obtained from the canopy-resolving simulations, 
$D$, to a form 
$D\approx C_{dh}U|U|$, where 
$U$ is the mean streamwise velocity. The drag coefficients obtained from cases T1, TP1 and T2 are portrayed in figure 3. This quadratic form provides a reasonable approximation of the drag force for 
$y^{+}\gtrsim 20$, once viscous effects are small. This is consistent with observations made in previous studies (Coceal et al. Reference Coceal, Thomas, Castro and Belcher2006, Reference Coceal, Thomas and Belcher2008; Böhm et al. Reference Böhm, Finnigan, Raupach and Hughes2013).
For the simulations labelled with the suffix ‘-H’, the presence of the canopy is replaced by a force 
$C_{dh}u_{i}|u_{i}|$ applied homogeneously below the canopy tips. This is the conventional homogeneous-drag model. It also requires the prescription of drag coefficients in the spanwise and wall-normal directions. We estimate these by rescaling the streamwise drag coefficient based on the relative change in the ‘blockage ratio’ of the canopy elements in the different directions (Luhar, Rominger & Nepf Reference Luhar, Rominger and Nepf2008), in the spirit of the method proposed by Luhar & Nepf (Reference Luhar and Nepf2013). The blockage ratio in each direction is proportional to the frontal area of the canopy elements in that direction. In the wall-normal and spanwise directions, this would be the top-view and the side-view areas, respectively. For the wall-normal drag, this assumption is particularly coarse, but Busse & Sandham (Reference Busse and Sandham2012) have shown that the flow is relatively insensitive to moderate changes in the wall-normal drag coefficient.
In the simulations labelled with the suffix ‘-H0’, a forcing 
$C_{dh}U|U|$ is applied in the region below the canopy tips, where 
$U(y)$ is the mean-velocity profile. The drag is only applied to the mean streamwise velocity, and has no fluctuating component. While the drag force in cases labelled ‘-H’ varies along any given wall-parallel plane depending on the local velocity, in cases labelled ‘-H0’ the drag force is homogeneous along any given wall-parallel plane, as it depends only on the mean velocity and the drag coefficient at that height. Note that, as the aforementioned drag models do not resolve the shape of the canopy elements, they also cannot capture the element-induced flow. In order to capture a part of the effect of the element-induced flow, the simulation labelled with the suffix ‘-L’ applies a drag 
$C_{dh}U|U|$, as in cases H0, but distributed in a reduced-order representation of the canopy elements. This representation consists of a 24-mode, 
$x{-}z$ Fourier truncation of the canopy geometry.
2.3 Effect of Reynolds number
Figure 4. The r.m.s. velocity fluctuations scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, from: red ——, case P2; blue ——, case 
$\text{P2I}_{Re}$; and dark green ——, case 
$\text{P2O}_{Re}$. In panels (a,c,e), the wall-normal coordinate is scaled in friction units, and in (b,d,f) in outer units. The black lines represent the smooth-wall simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$ and 1000. The smooth-wall data at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$ are taken from Lee & Moser (Reference Lee and Moser2015).
Figure 5. Viscous and Reynolds shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, from: red ——, case P2; blue ——, case 
$\text{P2I}_{Re}$; and dark green ——, case 
$\text{P2O}_{Re}$. In panel (a), the wall-normal coordinate is scaled in friction units, and in (b) in outer units. The black lines represent the smooth-wall simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$ and 1000. The smooth-wall data at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$ are taken from Lee & Moser (Reference Lee and Moser2015).
Figure 6. Drag coefficients, 
$C_{dh}=D/U^{2}$ for red ——, case P2; blue ——, case 
$\text{P2I}_{Re}$; and dark green ——, case 
$\text{P2O}_{Re}$.
In order to analyse the effect of the Reynolds number on the DNS results presented in the subsequent sections, we first compare the results of case P2 to those of cases 
$\text{P2I}_{Re}$ and 
$\text{P2O}_{Re}$. The simulations P2 and 
$\text{P2I}_{Re}$ have the same canopy parameters in friction units and 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$ and 1000, respectively. The root-mean-square (r.m.s.) velocity fluctuations and the Reynolds shear stresses for these simulations essentially collapse up to a height 
$y^{+}\approx 20$, as shown in figures 4(a,c,e) and 5(a). At heights larger than 
$y^{+}\gtrsim 20$, the magnitudes of velocity fluctuations and the Reynolds shear stresses are larger for case 
$\text{P2I}_{Re}$ than for case P2. This behaviour is consistent with that observed for smooth-wall flows at the corresponding Reynolds numbers (Sillero, Jiménez & Moser Reference Sillero, Jiménez and Moser2013), which suggests that the differences in the velocity fluctuations observed at heights 
$y^{+}\gtrsim 20$ do not result from the presence of the canopy. At heights of 
$y^{+}\gtrsim 200$, the velocity fluctuations of the canopy simulations collapse with those of their respective smooth-wall simulations, which indicates the recovery of outer-layer similarity. In addition, the effective canopy drag coefficients, 
$C_{dh}$, for the canopies of cases P2 and 
$\text{P2I}_{Re}$ also collapse in friction units, as shown in figure 6(a). These results show that the domain height used in case P2 is sufficiently large and does not constrain the flow within the canopy layer. Increasing the domain height further simply results in a larger region above the canopy layer exhibiting outer-layer similarity. This is consistent with the study of flows over cube canopies by Coceal et al. (Reference Coceal, Thomas, Castro and Belcher2006), who also noted that increasing the height of the domain beyond 
$\unicode[STIX]{x1D6FF}/\ell _{s}\approx 4$ did not have a significant effect on the flow within the canopy layer.
We now compare the results from cases P2 and 
$\text{P2O}_{Re}$, which have the same canopy parameters when scaled in outer units, but different friction Reynolds numbers. The canopy heights for both of these cases is 
$\ell _{s}/\unicode[STIX]{x1D6FF}\approx 0.2$, and in both cases the elements extend well into the logarithmic, self-similar region of the flow. Close to the wall, 
$y^{+}\lesssim 20$, the velocity fluctuations and Reynolds shear stresses for these cases collapse when scaled in friction units. At larger heights, 
$y^{+}>20$, the streamwise velocity fluctuations and the Reynolds shear stresses essentially collapse when the wall-normal coordinate is scaled in outer units, as shown in figures 4(b) and 5(b). The cross-velocity fluctuations in the region between 
$y^{+}\approx 20$ and 
$y/\unicode[STIX]{x1D6FF}\approx 0.2$ are larger for case 
$\text{P2O}_{Re}$ compared to case P2, as illustrated in figure 4(d,f). This increase, however, is also observed in the cross-velocity fluctuations of the corresponding smooth-wall simulations. Beyond a height of 
$y/\unicode[STIX]{x1D6FF}\approx 0.25$, we observe that the velocity fluctuation and Reynolds shear stress profiles for cases P2, 
$\text{P2O}_{Re}$, 
$\text{P2I}_{Re}$ and the smooth-wall simulations coincide. Furthermore, the effective canopy drag coefficient profiles for cases P2 and 
$\text{P2O}_{Re}$, portrayed in figure 6(b), also collapse when scaled in outer units, consistent with the observations of Cheng & Castro (Reference Cheng and Castro2002) for cube canopies. Therefore, the drag coefficient of the canopy is essentially independent of the Reynolds number, implying that the canopy is in the fully rough regime (Nikuradse Reference Nikuradse1933). These results therefore suggest that the conclusions drawn in the following sections from simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$ should also be relevant for higher-Reynolds-number flows. Note also that the quotient 
$D/U^{2}$ becomes close to constant for 
$y^{+}\gtrsim 25$, indicating that the total canopy drag, 
$D$, is essentially quadratic above this height.
3 Canopy-resolving simulations
In this section, we present and discuss the scaling of turbulent fluctuations in sparse canopies, and compare them with those over a smooth wall. Over a smooth wall, the balance of stresses within the channel can be obtained by averaging the momentum equations in the wall-parallel directions and time, and integrating in 
$y$, which yields
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{d}P}{\text{d}x}y+\unicode[STIX]{x1D70F}_{w}=-\overline{uv}+\unicode[STIX]{x1D708}\frac{\text{d}U}{\text{d}y}, & \displaystyle\end{eqnarray}$$ where 
$\unicode[STIX]{x1D70F}_{w}$ is the wall shear stress, 
$\text{d}P/\text{d}x$ is the mean streamwise pressure gradient, 
$-\overline{uv}$ is the Reynolds shear stress, 
$U$ is the mean streamwise velocity and 
$\unicode[STIX]{x1D708}$ is the kinematic viscosity. Particularising (3.1) at 
$y=\unicode[STIX]{x1D6FF}$, we obtain the expression for the wall shear stress and the friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$:
$$\begin{eqnarray}\displaystyle & \displaystyle u_{\unicode[STIX]{x1D70F}}^{2}=\unicode[STIX]{x1D70F}_{w}=-\unicode[STIX]{x1D6FF}\frac{\text{d}P}{\text{d}x}. & \displaystyle\end{eqnarray}$$In the presence of a canopy, the stress balance also includes the drag exerted by the canopy elements,
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{d}P}{\text{d}x}y+\unicode[STIX]{x1D70F}_{w}=-\overline{uv}+\unicode[STIX]{x1D708}\frac{\text{d}U}{\text{d}y}-\int _{0}^{y}D~\text{d}y, & \displaystyle\end{eqnarray}$$ where the canopy drag, 
$D$, averaged in 
$x$, 
$z$ and time is zero in the region above the canopy tips, 
$y>h$. This stress balance is typically used in the canopy literature to calculate the canopy drag stress (Dunn, Lopez & Garcia Reference Dunn, Lopez and Garcia1996; Ghisalberti & Nepf Reference Ghisalberti and Nepf2004). Equation (3.3) can be rewritten as
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{d}P}{\text{d}x}y+\unicode[STIX]{x1D70F}_{w}+\int _{0}^{h}D\,\text{d}y=-\overline{uv}+\unicode[STIX]{x1D708}\frac{\text{d}U}{\text{d}y}+\int _{y}^{h}D\,\text{d}y, & \displaystyle\end{eqnarray}$$ so that the net drag, 
$\unicode[STIX]{x1D70F}_{w}+\int _{0}^{h}D\,\text{d}y$, is on the left-hand side, as in (3.1). From this net drag, a ‘global’ friction velocity can be defined,
$$\begin{eqnarray}\displaystyle & \displaystyle u_{\unicode[STIX]{x1D70F}}^{2}=\unicode[STIX]{x1D70F}_{w}+\int _{0}^{h}D\,\text{d}y=-\unicode[STIX]{x1D6FF}\frac{\text{d}P}{\text{d}x}. & \displaystyle\end{eqnarray}$$ This is the equivalent of the smooth-wall 
$u_{\unicode[STIX]{x1D70F}}$ of (3.2) for canopy flows.
Figure 7. Schematic representation of the stress balance in a channel with canopy elements.
Figure 8. Stress profiles within the channel for case P1. (a) Dark green ——, full Reynolds shear stress; dark green – – –, background-turbulence Reynolds shear stress; red ——, drag stress; and blue ——, viscous stress; all scaled with 
$u_{\unicode[STIX]{x1D70F}}$. (b) Dark green ——, background-turbulence Reynolds shear stress; and blue ——, viscous stress; both scaled with 
$u^{\ast }$. The black lines represent the smooth-wall case, S.
While in smooth-wall flows the total stress is the sum of the viscous and Reynolds shear stresses alone and is linear in 
$y$, in canopy flows that linear sum has an additional contribution from the canopy drag as evidenced by the right-hand side of (3.4). This equation also portrays that, at any given height 
$y$, the sum of the streamwise shear stresses, 
$-\overline{uv}+\unicode[STIX]{x1D708}\,\text{d}U/\text{d}y$, and the drag from the canopy above that height, 
$\int _{y}^{h}D\,\text{d}y$, are balanced by the force exerted by the pressure gradient above. This can also be obtained from an integral balance of forces between heights 
$y$ and 
$\unicode[STIX]{x1D6FF}$, and is illustrated by the sketch in figure 7. Outside the canopy, the drag term is zero, and the magnitude of the viscous and Reynolds shear stresses is similar to that over smooth walls. Within the canopy, however, the canopy drag can dominate, and the viscous and Reynolds shear stresses are smaller than over smooth walls, as shown in figure 8(a).
Figure 9. The r.m.s. velocity fluctuations and shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, for: red ——, case TP1; and blue ——, case T1. Solid lines represent the full velocity fluctuations and dashed lines represent the background-turbulence fluctuations. The black lines represent the smooth-wall case, S, for reference.
As can be observed in figure 1, the presence of the canopy elements induces a coherent flow. Several studies have shown that the flow around the canopy elements and the flow far away from them have significantly different characteristics, and consequently they are typically studied separately (Finnigan Reference Finnigan2000; Bailey & Stoll Reference Bailey and Stoll2013; Böhm et al. Reference Böhm, Finnigan, Raupach and Hughes2013). A commonly used technique to separate the element-induced flow from the background turbulence is through triple decomposition (Reynolds & Hussain Reference Reynolds and Hussain1972)
$$\begin{eqnarray}\boldsymbol{u}=\boldsymbol{U}+\widetilde{\boldsymbol{u}}+\boldsymbol{u}^{\prime },\end{eqnarray}$$ where 
$\boldsymbol{u}$ is the full velocity. The mean velocity, 
$\boldsymbol{U}$, is obtained by averaging the flow in time and space. The element-induced velocity, also referred to as the dispersive flow, 
$\widetilde{\boldsymbol{u}}$, is obtained by ensemble-averaging the flow in time alone. We refer to the remaining part of the flow, 
$\boldsymbol{u}^{\prime }$, as the incoherent, background-turbulence velocity. Similarly, we refer to the Reynolds shear stress calculated using the full velocity, 
$\overline{uv}$, as the ‘full’ Reynolds shear stress and that calculated using the incoherent, background-turbulence velocity, 
$\overline{u^{\prime }v^{\prime }}$, as the background-turbulence Reynolds shear stress. The difference between these two stresses gives a measure of the element-induced, dispersive stress, 
$\overline{\tilde{u} \tilde{v}}=\overline{uv}-\overline{u^{\prime }v^{\prime }}$. Note that this is slightly different from the commonly used notation, where the dispersive and background-turbulence Reynolds shear stresses are treated distinctly and the ‘full’ Reynolds shear stress is not labelled (e.g. Coceal et al. Reference Coceal, Thomas, Castro and Belcher2006). Essentially identical to triple decomposition, ‘double averaging’ (Raupach & Shaw Reference Raupach and Shaw1982) can also be used to separate the element-induced and background-turbulence flows (Finnigan Reference Finnigan2000; Nepf Reference Nepf2012; Bai, Katz & Meneveau Reference Bai, Katz and Meneveau2015; Giometto et al. Reference Giometto, Christen, Meneveau, Fang, Krafczyk and Parlange2016; Yan et al. Reference Yan, Huang, Miao, Cui and Zhang2017).
The intensity of the element-induced flow can vary with the shape (Balachandar, Mittal & Najjar Reference Balachandar, Mittal and Najjar1997; Taylor et al. Reference Taylor, Palombi, Gurka and Kopp2011), permeability (Yu et al. Reference Yu, Zeng, Lee, Bai and Low2010; Ledda et al. Reference Ledda, Siconolfi, Viola, Gallaire and Camarri2018) and distribution of the canopy elements. It is possible, however, for canopies to have different element-induced flows but similar background turbulence. To illustrate this, we compare two canopy simulations, T1 and TP1, which have similar canopy layouts and roughly the same net drag, as shown in figure 9(d). The difference between the two cases is that in T1 the canopy elements are impermeable, whereas in TP1 some flow penetrates into the elements. The r.m.s. fluctuations of the full and background-turbulence velocity components for these cases are shown in figure 9. The magnitude of the full streamwise fluctuations within the canopy is significantly larger for T1 than for TP1. This increase, however, can be attributed essentially to the stronger element-induced fluctuations generated by the impermeable canopy elements. This is evidenced by the fact that the background-turbulence streamwise fluctuations for both cases essentially collapse, as shown in figure 9(a). The cross-flow fluctuations and Reynolds shear stress profiles for both these canopies are also similar. The impermeable canopy, however, has a slightly larger damping effect on the spanwise fluctuations.
Figure 10. The r.m.s. velocity fluctuations scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, in (a,c,e) and with the local friction velocity, 
$u^{\ast }$, in (b,d,f). The lines represent: red ——, case P0; violet ——, case P1; and blue ——, case P2. In panels (a,c,e), solid lines represent the full velocity fluctuations and dashed lines the background-turbulence fluctuations. In panels (b,d,f), only the background-turbulence fluctuations are portrayed. The black lines represent the smooth-wall case, S, for reference.
The fluctuating velocities portrayed in figure 10(a,c,e) are scaled using the ‘global’ friction velocity defined by (3.5), which includes the full contribution of the canopy drag. Tuerke & Jiménez (Reference Tuerke and Jiménez2013) studied smooth-wall flows with artificially forced mean profiles, and observed that the turbulent fluctuations in such flows scaled with the local sum of the viscous and Reynolds shear stresses, or the local stress 
$\unicode[STIX]{x1D70F}_{f}$, at each height. This was the case even when 
$\unicode[STIX]{x1D70F}_{f}$ was not linear with 
$y$ due to the artificial forcing. This idea has been expanded on by Lozano-Durán & Bae (Reference Lozano-Durán and Bae2019), who proposed that the energy-containing turbulent scales in the logarithmic layer in smooth-wall flows also scale with local velocity and length scales at each height, irrespective of the location of the wall. Tuerke & Jiménez (Reference Tuerke and Jiménez2013) defined a ‘local’ friction velocity, 
$u^{\ast }$, by linearly extrapolating the local stress at each height to the wall:
$$\begin{eqnarray}{u^{\ast }(y)}^{2}=\frac{\unicode[STIX]{x1D6FF}}{\unicode[STIX]{x1D6FF}-y}\unicode[STIX]{x1D70F}_{f}(y).\end{eqnarray}$$ Notice that, for a smooth unforced channel, 
$u^{\ast }=u_{\unicode[STIX]{x1D70F}}$ at every height. Following Tuerke & Jiménez (Reference Tuerke and Jiménez2013), we define the sum of the viscous and background-turbulence Reynolds stresses as the ‘fluid’ stress 
$\unicode[STIX]{x1D70F}_{f}$. In the present work, we only discuss the scaling of the background-turbulence fluctuations. Hence, only the contribution of the background-turbulence Reynolds shear stresses to 
$\unicode[STIX]{x1D70F}_{f}$ is considered. A similar concept was also proposed by Högström, Bergström & Alexandersson (Reference Högström, Bergström and Alexandersson1982) for flows over urban canopies. They scaled turbulence with a local friction velocity, defined as the square root of the magnitude of the local Reynolds shear stress, but had measurements only at heights where the contribution of the viscous stress to 
$\unicode[STIX]{x1D70F}_{f}$ would be small. Using 
$u^{\ast }$, a local viscous length scale can also be defined, 
$\unicode[STIX]{x1D708}/u^{\ast }$, and from it an effective viscous height, 
$y^{\ast }=yu^{\ast }/\unicode[STIX]{x1D708}$. Both 
$u^{\ast }$ and 
$y^{\ast }$ are portrayed, for the prismatic-post canopies, in figure 11. Near the canopy tips, where the element-induced drag is no longer present, the local friction velocity, 
$u^{\ast }$, becomes equal to the global 
$u_{\unicode[STIX]{x1D70F}}$, and 
$y^{\ast }$ becomes equal to 
$y^{+}$. Making the canopy sparser reduces the canopy drag, and hence the difference between 
$u^{\ast }$ and 
$u_{\unicode[STIX]{x1D70F}}$ within the canopy reduces with increasing canopy sparsity.
Figure 11. Variation of (a) 
$y^{\ast }$ and (b) 
$u^{\ast }$ with height for cuboidal post canopies, with sparsity increasing from red to blue: red ——, denser case P0; violet ——, intermediate case P1; blue ——, sparser case P2; and black ——, smooth-wall case, S.
Figure 12. Shear stresses scaled using (a) the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, and (b) the local friction velocity, 
$u^{\ast }$, for: red ——, case P0; violet ——, case P1; and blue ——, case P2. The solid and dashed lines in panel (a) represent the full and background-turbulence Reynolds shear stresses, respectively. In panel (b) only the background-turbulence Reynolds shear stresses are portrayed. The black lines represent the smooth-wall case, S.
When scaled with 
$u_{\unicode[STIX]{x1D70F}}$, as is done conventionally, the viscous and Reynolds shear stresses near the base of the canopy are highly damped compared to smooth walls. However, the balance of these stresses in 
$\unicode[STIX]{x1D70F}_{f}$ remains close to that over smooth walls. This is illustrated in figures 8(b) and 12(b), which portray the terms in the stress balance within a channel with canopies scaled with 
$u_{\unicode[STIX]{x1D70F}}$ and 
$u^{\ast }$. The similarity of the viscous and Reynolds shear stresses in the canopy-flow and smooth-wall cases suggests that the canopy acts on the background turbulence essentially through changing their local scale, rather than through a direct interaction of the canopy elements with the flow. To explore the scaling further, the background-turbulence r.m.s. fluctuations for the prismatic post canopies are portrayed scaled with 
$u^{\ast }$ in figure 10(b,d,f). Scaling the fluctuations with the conventional 
$u_{\unicode[STIX]{x1D70F}}$ shows a reduction of the fluctuations within the canopy compared to a smooth wall, as shown in figure 10(a,c,e). With our proposed scaling with 
$u^{\ast }$, in contrast, the streamwise fluctuations appear similar to those in a smooth channel. The increase in spanwise and wall-normal fluctuations, shown in figure 10(b,d,f), suggests, however, that there is a relative increase in the intensity of the cross-flow within the canopy compared to a smooth channel. The velocity fluctuations and the shear stresses for the higher-Reynolds-number simulations and for the T-shaped canopies also exhibit similar behaviour, and are provided in the Appendix for reference.
Although 
${u^{\prime }}^{\ast }$ and 
$\overline{u^{\prime }v^{\prime }}^{\ast }$ within the canopy appear similar to those over smooth walls, there are some differences in the distribution of energy across different length scales in the flow, particularly in the region close to the wall. In order to examine this, we compare the spectral energy densities, at 
$y^{\ast }=15$, for a smooth wall and for case P1 in figure 13. This is the height roughly corresponding to the location at which the magnitude of the fluctuations peaks in smooth-wall flows (Jiménez & Pinelli Reference Jiménez and Pinelli1999). In global units, the energy is observed to be in larger wavelengths when compared to a smooth channel, especially in 
$\unicode[STIX]{x1D706}_{z}$. In local scaling, however, there is a greater overlap of the regions with highest intensity, particularly for 
$E_{uu}$ and 
$E_{uv}$. In addition, the canopy case exhibits a concentration of energy at the canopy wavelengths and its harmonics. Note that the canopy spacing, for case P1, at 
$y^{\ast }\approx 15$ is reduced to 
$L_{x}^{\ast }=L_{z}^{\ast }\approx 100$, while in global scaling it is 
$L_{x}^{+}=L_{z}^{+}\approx 200$. The increase in the energy in the canopy wavelengths is a reflection of the element-induced flow. The large streamwise scales, in turn, are damped by the presence of the canopy, which results in a reduction of their energy.
Figure 13. Pre-multiplied spectral energy densities for case P1 (filled contours) and for case S (line contours) normalised with their respective r.m.s. values at (a–d) 
$y^{+}=15$ and (e–h) 
$y^{\ast }=15$. The contours, from the left to right columns, are in increments of 0.03, 0.06, 0.05 and 0.06, respectively.
Figure 14. Pre-multiplied spectral energy densities for case P1 (filled contours) and for case S (line contours) at 
$y^{\ast }=105$, normalised by their respective 
$u^{\ast }$. The contours in (a–d) are in increments of 0.125, 0.06, 0.075 and 0.06, respectively.
Figure 15. Pre-multiplied spectral energy densities at 
$y^{+}=250$: red ——, case P0; violet ——, case P1; and blue ——, case P2. All are normalised by their respective 
$u_{\unicode[STIX]{x1D70F}}$. Filled contours represent case S. The contours in (a–d) are in increments of 0.075, 0.04, 0.06 and 0.03, respectively.
Figure 16. Mean-velocity profiles, from the canopy-resolving simulations. Lines represent: red ——, case T1; violet ——, case TP1; blue ——, case T2; red – – –, case P0; violet – – –, case P1; and blue – – –, case P2. Here 
$U_{c}$ is the mean velocity at 
$y=\unicode[STIX]{x1D6FF}$. The black lines represent the smooth-wall case, S.
Figure 17. Pre-multiplied spectral energy densities at 
$y^{\ast }\approx 115$, normalised by their respective r.m.s. values. The line contours represent: (a–d) case P2; (e–h) case P1; (i–l) case P0; and (m–p) case PD. Filled contours represent case S. The contour increments in the leftmost to rightmost columns are 0.029, 0.048, 0.045 and 0.045, respectively.
The differences in the energy distribution observed within the canopy eventually disappear above it. To illustrate this, figure 14 portrays the spectra near the canopy tips, 
$y^{\ast }\approx 105$. Here, the concentration of energy in the canopy wavelengths and its harmonics is weak, and the smaller scales in the flow are smooth-wall-like. There is, however, still a deficit of energy in large streamwise wavelengths compared to a smooth wall, associated with the damping of these scales by the canopy elements, as discussed in the previous paragraph. This effect diminishes away from the canopy, and the spectra are essentially smooth-wall-like for 
$y^{+}\gtrsim 250$, as shown in figure 15, indicating that outer-layer similarity is recovered beyond this height. Consequently, we can also conclude that this height marks the extent of the roughness sublayer of the canopies. The recovery of outer-layer similarity is also reflected in the mean-velocity profiles of the canopy-flow simulations, portrayed in figure 16, which exhibit logarithmic-law behaviour with a standard Kármán constant when shifted by a suitable displacement height, 
$d$ (Jackson Reference Jackson1981).
So far, we have mainly focused on the results for case P1, with prismatic canopy elements with spacings 
$L_{x}^{+}=L_{z}^{+}\approx 200$. We now discuss the effect of the canopy element geometry and spacing. An increased sparsity results in an increase in the magnitude of both the full and background-turbulence velocity fluctuations, as shown in figure 10. In local scaling, however, the background-turbulence fluctuations follow a similar trend to that observed for case P1. We observe that 
${u^{\prime }}^{\ast }$ and 
$\overline{u^{\prime }v^{\prime }}^{\ast }$ appear smooth-wall-like, while there is a relative increase in the magnitude of the cross-fluctuations compared to those over a smooth wall. For the denser canopy of case P0, on the other hand, the fluctuations become less similar to those over smooth walls. The streamwise fluctuations are damped more intensely within the canopy, and there are additional Reynolds shear stresses near the wall. Figure 17 shows that, compared to the sparser canopies, P0 has an accumulation of energy in streamwise wavelengths corresponding to the canopy harmonics but across a range of spanwise wavelengths. These regions of excess energy have also been noted by Abderrahaman-Elena, Fairhall & García-Mayoral (Reference Abderrahaman-Elena, Fairhall and García-Mayoral2019), who studied densely packed cuboidal roughness. They noted that these regions were an imprint of the large, background-turbulence scales modulating the smaller-scale coherent flow generated by the roughness. This effect diminishes as the canopy element spacing is made larger than the energetic scales in the background turbulence, as evidenced by the lack of these regions in the spectra of the sparser canopies portrayed in figures 13, 14 and 17. For case P2, the spectra are already close to smooth-wall-like near the canopy tips, suggesting that both the element-induced flow and the damping of large scales are already weak at this height.
The results discussed above suggest that there is a progressive departure from smooth-wall-like behaviour in canopy flows as the element spacing is reduced, which is consistent with the observations of Poggi et al. (Reference Poggi, Porporato, Ridolfi, Albertson and Katul2004) and Huang et al. (Reference Huang, Cassiani and Albertson2009). If the element spacing was reduced even further, eventually we would expect the complete breakdown of smooth-wall-like dynamics within the canopy. In the resulting dense canopy, the flow near the canopy tips would be characterised by the Kelvin–Helmholtz-like, mixing-layer instability (Raupach et al. Reference Raupach, Finnigan and Brunet1996; Finnigan Reference Finnigan2000; Poggi et al. Reference Poggi, Porporato, Ridolfi, Albertson and Katul2004; Nepf Reference Nepf2012). To investigate this effect, we have conducted an additional simulation of the prismatic post canopies in an even denser arrangement, case PD. In order to contrast the flow characteristics of sparse canopies with those of dense ones, we now compare the results from the simulation PD to that of the sparsest canopy studied, case P2. Deep within the dense canopy of case PD, the viscous and the Reynolds shear stresses are negligible and the drag stress dominates, as shown in figure 18(d). While for the sparse canopy of case P2, the magnitude of the velocity fluctuations near the canopy tips is similar to that over a smooth wall, for the dense canopy of case PD their magnitude is considerably reduced. Furthermore, the magnitude of the element-induced streamwise fluctuations in the dense canopy is negligible compared to that of the sparse canopy, where it is observed to constitute up to 30 % of the total fluctuations, similar to the observations of Poggi & Katul (Reference Poggi and Katul2008).
Figure 18. The r.m.s. velocity fluctuations and shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, from: red ——, case PD; and blue ——, case P2. Solid lines represent the full velocity fluctuations and dashed lines represent the background-turbulence fluctuations. The vertical dashed line marks the location of the canopy tip plane. The black lines represent the smooth-wall case, S, for reference.
Figure 19. Instantaneous realisations of the wall-normal velocity at 
$y^{+}\approx 120$, normalised by 
$u_{\unicode[STIX]{x1D70F}}$. Panels (a)–(e) represent cases S, P2, P1, P0 and PD, respectively. The clearest and darkest colours in (a)–(e) represent intensities of 
$\pm (1.5,1.5,1.5,1,1)$, respectively.
Previous studies have noted the formation of Kelvin–Helmholtz-like instabilities near the canopy tips in dense canopies (Finnigan Reference Finnigan2000; Poggi et al. Reference Poggi, Porporato, Ridolfi, Albertson and Katul2004; Nepf Reference Nepf2012). When present, these instabilities leave a distinct footprint in 
$E_{vv}$ and 
$E_{uv}$, causing an increase in energy in a narrow range of streamwise wavelengths and for large spanwise wavelengths (García-Mayoral & Jiménez Reference García-Mayoral and Jiménez2011; Gómez-de Segura, Sharma & García-Mayoral Reference Gómez-de Segura, Sharma and García-Mayoral2018; Abderrahaman-Elena et al. Reference Abderrahaman-Elena, Fairhall and García-Mayoral2019). For the canopies of cases P0, P1 and P2, such a footprint is not observed in the spectral energy densities, portrayed in figure 17. In the energy densities of the wall-normal velocity for the dense canopy, case PD, we observe some concentration of energy in a range of streamwise wavelengths, 
$\unicode[STIX]{x1D706}_{x}^{+}\in 200{-}500$, consistent with the presence of a Kelvin–Helmholtz-like instability. The gradual breakdown of the smooth-wall-like behaviour of flows over canopies with decreasing element spacing can also be observed in instantaneous realisations of the wall-normal velocity, portrayed in figure 19. Dense canopies are not the focus of the present study, so the results from case PD are not discussed extensively here. Increasing the canopy density yet further can result in a more distinct imprint of the Kelvin–Helmholtz-like instability near the canopy tips, as discussed in Poggi et al. (Reference Poggi, Porporato, Ridolfi, Albertson and Katul2004) and Sharma & García-Mayoral (Reference Sharma and García-Mayoral2019). Further details about the formation and development of this instability over dense canopies can be found in Raupach et al. (Reference Raupach, Finnigan and Brunet1996), Ghisalberti & Nepf (Reference Ghisalberti and Nepf2002), Finnigan et al. (Reference Finnigan, Shaw and Patton2009) and Bailey & Stoll (Reference Bailey and Stoll2016).
4 Simulations with artificial forcing
The results discussed so far suggest that sparse canopies affect their surrounding flow through two mechanisms: an element-induced flow, and a change in the local scale for the background-turbulence fluctuations. With respect to the second mechanism, the effect of the canopy elements would be indirect, through modifying the mean-velocity profile and thus the local stress, 
$\unicode[STIX]{x1D70F}_{f}$. The latter would, in turn, set the scale for the fluctuations. If this is the case, applying the mean drag produced by the canopy on the mean flow alone should capture the essential effects of the canopy on the background turbulence. We test this in the simulations labelled with the suffix ‘-H0’. For the cases P0 and T1/TP1, we also compare the mean-only-drag simulations with conventional, homogeneous-drag simulations labelled with the suffix ‘-H’. These could be expected to be better representations when the canopy is dense enough for all the turbulent scales to perceive it in a homogenised fashion. Note that simulations T1-H0 and T1-H correspond to both the permeable and impermeable canopies of T1 and TP1, as they have similar net drags and drag coefficients.
Figure 20. Background-turbulence r.m.s. velocity fluctuations and shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, of the canopy-resolving and mean-only/homogeneous-drag simulations. The lines represent: red ——, case P0; red – – –, case P0-H0; red 
$\cdots \cdots$, case P0-H; violet ——, case P1; violet – – –, case P1-H0; blue ——, case P2; and blue – – –, case P2-H0. The black lines represent the smooth-wall case, S.
Figure 21. Background-turbulence r.m.s. velocity fluctuations and shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, of the canopy-resolving and mean-only/homogeneous-drag simulations. The lines represent: red ——, case T1; red – – –, case T1-H0; red 
$\cdots \cdots$, case T1-H; violet ——, case TP1; blue ——, case T2; and blue – – –, case T2-H0. The black lines represent the smooth-wall case, S.
The streamwise fluctuations and the Reynolds shear stresses of the mean-only-drag simulations are in good agreement with the corresponding background-turbulence fluctuations from the resolved canopies, except for case P0, as shown in figures 20 and 21. For the sparsest canopies, cases P2 and T2, the cross-flow fluctuations, particularly in the wall-normal direction, are slightly larger than their mean-only-drag counterparts. A likely reason for this discrepancy is the presence of an unsteady element-induced flow for these canopies, whose contribution cannot be filtered out by the conventional triple-decomposition technique that we have used. The fluctuating velocities scaled by the local friction velocities for the drag-force simulations are provided in the Appendix for reference.
For case P0, a homogeneous drag provides a better representation of the cross-flow fluctuations than the mean-only drag, and the streamwise fluctuations are not well represented by either forcing method. For this case, there is significant interaction between the element-induced flow and the background turbulence, as discussed in § 3. Thus, it is not surprising that neither the mean-only drag nor the homogeneous drag is able to capture the full effect of this canopy on the background turbulence.
Figure 22. Pre-multiplied spectral energy densities at (a–d) 
$y^{\ast }=15$ and (e–h) 
$y^{\ast }=105$, normalised by their respective 
$u^{\ast }$. Filled contours represent case TP1. Lines represent: blue ——, case T1-H0; red ——, case T1-H; and violet ——, case TP1-L. The contours in panels (a–h) are in increments of 0.3, 0.075, 0.175, 0.075, 0.12, 0.06, 0.075 and 0.05, respectively.
For the sparser canopies of T1 and TP1, compared to a mean-only drag, the homogeneous drag tends to overdamp the fluctuations within the canopy, particularly in the streamwise direction, as can be observed in figure 21(a). The excessive damping of fluctuations by a homogeneous drag, in comparison to a resolved canopy, was also noted by Yue et al. (Reference Yue, Parlange, Meneveau, Zhu, van Hout and Katz2007) and Bailey & Stoll (Reference Bailey and Stoll2013). Figure 22 shows that this decrease in the intensity of fluctuations within the canopy is mainly a result of damping of the smaller streamwise wavelengths in the flow, 
$\unicode[STIX]{x1D706}_{x}^{\ast }\lesssim 200$. The homogeneous-drag simulation, T1-H, reproduces well the larger scales of the resolved canopy simulation, TP1. This suggests that scales much larger than the canopy spacing still perceive the canopy as homogeneous. Using the mean-only drag recovers some of the smaller streamwise scales, but it does not act directly on the larger scales as the actual canopy does. As the element spacing is increased, the range of scales affected in a homogenised fashion is shifted to larger scales, so that the energetic turbulent scales become less damped. This is consistent with the results portrayed in figure 17, which show that, near the canopy tips, the dense canopy P0 damps the energy at 
$\unicode[STIX]{x1D706}_{x}^{\ast }\sim 1000{-}2000$, compared to smooth walls, while the sparse canopy P2 leaves these scales relatively undisturbed.
Figure 23. Drag force distribution in the streamwise direction in a plane passing through the canopy heads for case TP1-L (blue ——). Distribution of the mean-only drag force, as in case T1-H0 (red – – –). The location of the canopy elements is sketched in black at the bottom of the figure.
The accumulation of energy in the length scales of the order of the canopy wavelengths and its harmonics observed in the canopy-resolving simulations requires a discrete representation of the canopy elements. Hence, it cannot be captured by either the mean-only-drag or the homogeneous-drag approaches. To introduce information of the canopy layout in the model, we distribute the drag calculated from the mean flow into a reduced-order representation of the canopy elements in case TP1-L. The representation consists of a truncation in Fourier space in 
$x$ and 
$z$, of the actual layout. The procedure is illustrated in figure 23 by the streamwise distribution of the drag force used by this model. In addition to capturing the local scaling of the flow, discussed in § 4, this model is also able to represent the concentration of energy in the canopy scales for case TP1, as observed in figure 22. This is reflected by the collapse of the r.m.s. fluctuations of the full velocity components of TP1-L and TP1, as shown in figure 24. The magnitude of spanwise velocity fluctuations within the canopy of TP1-L is slightly larger than TP1, likely due to the fact that TP1-L does not apply any form of spanwise drag force. The drag force, although only applied in the streamwise direction, is also able to reproduce the canopy harmonics in the spectra of 
$E_{vv}$ and 
$E_{ww}$, which are caused by the deflection of the streamwise flow around the canopy elements as a result of continuity. The large scales in the flow are similar to those in the mean-only drag, as the drag in this case also does not act on these scales directly. This method, however, is only able to capture the weak coherent flow generated by the permeable canopy, and still underpredicts the full streamwise velocity fluctuations of the impermeable canopy.
Figure 24. The r.m.s. fluctuations and shear stresses of the full flow, scaled with the global friction velocity 
$u_{\unicode[STIX]{x1D70F}}$: blue ——, case TP1; red ——, case T1; red – – –, case TP1-H0; and violet ——, case TP1-L. The black lines represent case S.
5 Conclusions
In the present work, we have studied turbulent flows within and over sparse canopies. Two different canopy element geometries have been studied, each for various different element spacings. We have also compared canopies with permeable and impermeable elements in the same arrangement. The effect of the Reynolds number has also been examined by comparing the results from simulations of the same canopy in both inner and outer scalings at 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$ and 1000. The flow was decomposed into an element-induced component and a background-turbulence component. It was found that, although the element-induced flow in the permeable and impermeable canopies studied here differ, the background turbulence was essentially the same. A new scaling for the background-turbulence fluctuations within sparse canopies was proposed. This scaling uses the friction velocity based on the local sum at each height of the viscous and Reynolds shear stresses, 
$\unicode[STIX]{x1D70F}_{f}$, rather than the conventional friction velocity, based on the net drag. When scaled with the proposed local friction velocity, the background-turbulence fluctuations and the viscous and Reynolds shear stresses appear more smooth-wall-like, compared to when conventional total-drag scaling is used. This suggests that the sparse canopy acts in a large part on the background turbulence through a change in the local scale, rather than through a direct interaction with the canopy elements. Based on the proposed scaling, we investigated the extent to which a drag force acting only on the mean flow captures the effect of the canopy on the background turbulence. The mean-only drag directly modifies the mean flow alone, which in turn sets 
$\unicode[STIX]{x1D70F}_{f}$ and, hence, the scale for the fluctuations. We show that the mean-only drag is able to capture the background-turbulence fluctuations within the canopies better than a conventional, homogeneous drag. Neither approach is, however, sufficient to capture the element-induced flow. The latter can be partially recovered by redistributing the mean-only drag in a low-order representation of the canopy.
Acknowledgements
A.S. was supported by an award from the Cambridge Commonwealth, European and International Trust. Computational resources were provided by the ‘Cambridge Service for Data Driven Discovery’ operated by the University of Cambridge Research Computing Service and funded by EPSRC Tier-2 grant EP/P020259/1. Additional computing time was also provided by PRACE DECI-15. This work was also partly supported by the European Research Council through the Third Multiflow Summer Workshop. We would also like to thank Professor H. Nepf for her helpful feedback on an earlier version of this manuscript.
Declaration of interests
The authors report no conflicts of interest.
Appendix. Turbulence statistics in local scaling
The turbulent velocity fluctuations and Reynolds shear stresses for the simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$, 
$\text{P2I}_{Re}$ and 
$\text{P2O}_{Re}$, are compared with those from case P2 in both global and local scaling in figure 25.
Figure 25. The r.m.s. velocity fluctuations and shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, in (a,c,e,g) and with the local friction velocity, 
$u^{\ast }$, in (b,d,f,h). The lines represent: red ——, case P2; blue ——, case 
$\text{P2I}_{Re}$; and dark green ——, case 
$\text{P2O}_{Re}$. In (a,c,e,g), solid lines represent the full velocity fluctuations and dashed lines represent the background-turbulence fluctuations. In (b,d,f,h), only the background-turbulence fluctuations are portrayed. The solid black lines represent the smooth-wall simulations at 
$Re_{\unicode[STIX]{x1D70F}}\approx 520$ and 1000. The smooth-wall data at 
$Re_{\unicode[STIX]{x1D70F}}\approx 1000$ are taken from Lee & Moser (Reference Lee and Moser2015).
Figure 26. The r.m.s. velocity fluctuations and shear stresses scaled with the global friction velocity, 
$u_{\unicode[STIX]{x1D70F}}$, in (a,c,e,g) and with the local friction velocity, 
$u^{\ast }$, in (b,d,f,h). The lines represent: red ——, case T1; violet ——, case TP1; and blue ——, case T2. In (a,c,e,g), solid lines represent the full velocity fluctuations and dashed lines represent the background-turbulence fluctuations. In (b,d,f,h), only the background-turbulence fluctuations are portrayed. The black lines represent the smooth-wall case, S, for reference.
Figure 27. Background-turbulence r.m.s. velocity fluctuations and shear stresses, scaled with the local friction velocity, 
$u^{\ast }$, of the resolved canopy and mean-only/homogeneous-drag simulations. In the left-hand column, lines represent: red ——, case P0; red – – –, case P0-H0; red 
$\cdots \cdots$, case P0-H; violet ——, case P1; violet – – –, case P1-H0; blue ——, case P2; and blue – – –, case P2-H0. In the right-hand column, lines represent: red ——, case T1; red – – –, case T1-H0; red 
$\cdots \cdots$, case T1-H; violet ——, case TP1; blue ——, case T2; and blue – – –, case T2-H0. The black lines represent case S.
The turbulence statistics, in both global and local scaling, for the T-shaped canopies are portrayed in figure 26.
Figure 27 compares the background-turbulence velocity fluctuations from the canopy-resolving simulations with their corresponding mean-only-drag or homogeneous-drag simulations, in the proposed local scaling.
