the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
A onedimensional urban flow model with an eddydiffusivity massflux (EDMF) scheme and refined turbulent transport (MLUCM v3.0)
Negin Nazarian
Melissa Anne Hart
E. Scott Krayenhoff
Alberto Martilli
In recent years, urban canopy models (UCMs) have been used as fully coupled components of mesoscale atmospheric models as well as offline tools to estimate temperature and surface fluxes using atmospheric forcings. Examples include multilayer urban canopy models (MLUCMs), where the vertical variability of turbulent fluxes is calculated by solving prognostic momentum and turbulent kinetic energy (TKE, k) using mixing length scale (l) and drag parameterizations. These parameterizations are based on the wellestablished 1.5order k−l turbulence closure theory and are often informed by microscale fluid dynamics simulations. However, this approach can include simplifications such as assuming the same diffusion coefficient for momentum, TKE, and scalars. In addition, the dispersive stresses arising from spatially averaged flow properties have been parameterized together with the turbulent fluxes despite being controlled by different mechanisms. Both of these assumptions impact the quantification of the turbulent exchange of flow properties and subsequent air temperature predictions in urban canopies. To assess these assumptions and improve corresponding parameterization, we analyzed 49 largeeddy simulations (LES) for idealized urban arrays, encompassing variable building height distributions and a comprehensive range of urban densities (${\mathit{\lambda}}_{p}\in [\mathrm{0.0625},\mathrm{0.64}]$) seen in global cities. We find that the efficiency of turbulent transport (numerically described via diffusion coefficients) is similar for scalars and momentum but is 3.5 times higher for TKE. Additionally, parameterizing the dispersive momentum flux using the k−l closure was a source of error, while scaling with the pressure gradient and urban morphological parameters appears more appropriate. In response to these findings, we propose two changes to the previous version of MLUCM: (a) separate characterization for turbulent diffusion coefficient for momentum and TKE and (b) introduction of an explicit physicsbased “massflux” term to represent the fraction of the dispersive momentum transport directly induced from buildings as an amendment to the existing “eddydiffusivity” framework. The updated onedimensional model, after being tuned for building height variability, is further compared against the original LES results and demonstrates improved performance in predicting vertical turbulent exchange in urban canopies.
 Article
(9882 KB)  Fulltext XML

Supplement
(2703 KB)  BibTeX
 EndNote
The urban canopy is a unique and complex land cover type in climate models (Oke et al., 2017). As climate models become more capable of highresolution simulations, the challenge resides in developing and validating urban canopy models (UCMs) to account for finerscale interactions with complex urban forms. Applications such as running UCMs to accurately describe urban surface fluxes to atmospheric boundary layers (Schoetter et al., 2020) and UCMs as offline models coupling with mesoscale forcing to simulate incanopy processes (Redon et al., 2020) can then benefit from such infrastructure developments.
The development of UCMs relies on the similarity between flow through the urban canopy layer (UCL) and other more commonly investigated types of flow. For example, Macdonald (2000) modified the exponential wind profile from flow over vegetative canopy (Finnigan, 2000) based on the aerodynamic roughness properties of the urban surface to predict the urban wind speed that is commonly used in singlelayer urban canopy models (SLUCM (Kusaka et al., 2001)). The 1.5order turbulence closure scheme (Bougeault and Lacarrere, 1989) used for turbulence parameterization in the planetary boundary layer (PBL) was also extended to the UCL, facilitating the development of the multilayer urban canopy model (MLUCM (Nazarian et al., 2020)). The urban form and fabric, however, feature distinct behaviors that require modifications to length scales (Cheng and PortéAgel, 2021) and form drag (Sützl et al., 2021b) that are uniquely developed to relate urban morphological parameters (Lu et al., 2023b). In addition, some physical processes commonly seen in urban environments, such as thermal instability, also result in unique flow characteristics in urban canopies that have been separately evaluated in urban canopy models (Santiago et al., 2014; Simónmoral et al., 2016; Nazarian and Kleissl, 2015).
The modeling of urban surfaces develops along with the capability of computational power and availability of realistic building morphology (Sirko et al., 2021). As the community is pushing the resolution of regional climate models to the subkilometer scale (Sützl et al., 2021a), there is an urgent need for a more detailed characterization of finer(subgrid)scale (and higherorder) urban phenomena (such as building wake and canyon recirculation (Oke et al., 2017)). The most common way to enable UCMs to capture these factors is by relating the variation of urban geometric parameters to the corresponding flow statistics based on computational fluid dynamics (CFD) models, where the largeeddy simulation (LES) offers a feasible way to obtain the higher moment of flow statistics.
The development of UCMs includes designing modules to predict flow properties such as wind speed (Castro, 2017), turbulent kinetic energy (TKE) (Christen et al., 2009), and scalars (Lim et al., 2022) that interact with urban structures differently. For example, increasing the frontal density of roughness elements was found to amplify the turbulent momentum transport but hinder the turbulent scalar transport (Li and BouZeid, 2019). However, efforts to identify and characterize the dissimilarity of the transport mechanisms of flow properties are rare due to the lack of experimental evidence to support quantitative characterizations. Therefore, momentum, TKE, and scalar diffusion efficiency in modeling practice are typically assumed to be the same (Martilli et al., 2015; Nazarian et al., 2020), introducing errors that may degrade model performance.
Multilayer UCMs are fundamentally more versatile than singlelayer UCMs in characterizing the urban effects (Hendricks et al., 2020; Lipson et al., 2023) and incanopy processes, and the use of UCMs has evolved into different implementations (Kondo et al., 2005; Yuan et al., 2019; Santiago and Martilli, 2010b). The building effects parameterization (BEP (Martilli et al., 2002)) model and developments to its flow parameterization through the multilayer urban canopy models (MLUCMs (Santiago and Martilli, 2010a)) have been directed toward multiple applications such as simulation of building–tree interactions (BEPtree (Krayenhoff et al., 2020)) and indoor building energy exchanges with outdoor climate (BEPBEM (Salamanca et al., 2010)). Simulations with MLUCMs are particularly sensitive to details of the parameterization that influence the strength of turbulence intensity within the canopy that facilitates turbulent transport of flow properties. Thus, accurate quantification of TKE is essential; however, previous work leading to the development of MLUCM v2.0 has shown an underestimation of canopy TKE (Nazarian et al., 2020). Moreover, canopy TKE must be accurately quantified for simplified neighborhoods before the addition of the flow effects arising from the added complexity of natural urban areas, such as subfacetscale structures (overhangs, awnings, irregular building shapes, etc.) or trees (Krayenhoff et al., 2020).
Arising from spatially averaged flow properties (that is core to MLUCM development), dispersive fluxes illustrate the transport of variables by timemean structures smaller than the averaging grid size, constituting another unique urban phenomenon (Poggi and Katul, 2008b). The strength of the dispersive flux is highly related to the horizontal (e.g., Lu et al., 2023b) and vertical (e.g., Xie et al., 2008) urban structures and exhibits substantial spatial variability (Harman et al., 2016) that precludes generalized characterizations. Therefore, dispersive flux has only been pragmatically considered in Simónmoral et al. (2016) based on the Ktheory framework as an increment to its turbulent counterpart. In practice, this approach resulted in an enhanced eddy diffusivity and implicitly assumed the correlation between dispersive flux and turbulence intensity and the vertical gradient of the flow property. However, the physical basis of such a correlation is still unclear and may not comprehensively improve the predictability of UCMs. Varying with the shape of roughness elements, the dispersive component can be about 50 % of the total transport over noncubic blocks (Li and BouZeid, 2019). Furthermore, realistic urban neighborhoods induce even higher dispersive contributions due to the more substantial spatial flow variability (Akinlabi et al., 2022; Giometto et al., 2016). Consequently, the lack of dispersive stress characterization restricts the study of canopy flows to relatively less complex and homogeneous flow, such as on gentle topography (e.g., Finnigan and Belcher, 2004) with slowly varying packing density (Ross, 2012).
Based on the aforementioned fundamental deficiencies embedded in the multilayer model, a better performance of UCMs on more complex flows requires a finer and more physically based characterization of dispersive fluxes. In this research, based on a recently developed urban flow dataset (Lu et al., 2023a), we aim to assess and improve the prediction of the TKE and momentum in the urban canopy through a more accurate parameterization of turbulent transport in the MLUCM (Nazarian et al., 2020) in two ways:

separate characterization of transport efficiency of different flow properties.

introduce a physicsbased modeling of the dispersive transport of momentum.
This paper first describes the flow dataset using the Parallelized LargeEddy Simulation Model (PALM) in Sect. 2.1. Then different components of MLUCM v2.0 (Nazarian et al., 2020) are revisited for the conventional eddy diffusion module in Sect. 3.1 and for the drag and dissipation module in Sect. 3.3. A novel physically based massflux parameterization is proposed based on the pressure gradient scaling and the geometric condition of the neighborhood in Sect. 3.2. In Sect. 4, the proposed changes for MLUCM are tested and compared across 49 urban arrays with different horizontal (staggered, aligned) and vertical (height variabilities) arrangements. Section 5 summarizes the findings of this study and provides perspectives for future developments of the multilayer model.
In this study, we used two methods to assess the turbulent transport behavior of the urban canopy flow, i.e., computational fluid dynamics (CFD) simulations and an offline multilayer urban canopy model. In Sect. 2.1, with the help of LES simulations recording the thirdorder moments in turbulent flows, we evaluate differences in turbulent transport of momentum, TKE, and scalars. Following the differences, we further revisit the turbulent parameterization in Sects. 3.1 and 3.2 followed by a comprehensive evaluation of model performance against LES in Sect. 4.
2.1 Largeeddy simulation (LES) setup and dataset
Conventionally, two configurations of idealized building arrays, “aligned” and “staggered” (Coceal et al., 2006), are employed to approximate horizontal arrangements of real neighborhoods. The aligned array simulates the situation in natural communities where streets align with the prevailing wind direction. The staggered array, by contrast, is mainly employed to calibrate UCM (Santiago and Martilli, 2010b; McNorton et al., 2021) in that it has no significant corridors (i.e., streets) aligned with the wind, potentially resembling the frequent situation where wind and street directions do not align. Both choices simplify realistic urban neighborhoods, but staggered arrays provide a closer approximation to average conditions in real cities (Nazarian et al., 2020). Accordingly, the staggered configuration was selected to obtain urban wind flow in determining parameterizations in MLUCM.
Simulations are conducted in the Parallelized LargeEddy Simulation Model (PALM, version r4554) (Maronga et al., 2020) with λ_{p} ranging from 0.0625 to 0.64, representing urban densities seen in global cities. The numerical setup follows (Lu et al., 2023a; Nazarian et al., 2020; Lu et al., 2023b), which has been validated against direct numerical simulation (DNS) (Coceal et al., 2007) and wind tunnel experiments (Brown et al., 2001) but will be listed here for completeness.
The computational domain is discretized equidistantly using secondorder central differences (Piacsek and Williams, 1970) with a horizontal and staggered Arakawa Cgrid in the vertical direction. The momentum field is solved using the filtered Navier–Stokes equation with the Boussinesq approximation with the time integration following a minimal storage scheme (Williamson, 1980). The covariance terms from the filtering procedure were parameterized using a 1.5order closure scheme after (Deardorff, 1980). The pressure perturbation was calculated in Poisson's equation and was solved by the FFTW scheme (Frigo and Johnson, 1998). Apart from solving momentum and pressure fields, a passive scalar was introduced with an initial surface value c_{0}=0.06 and a constant negative gradient of $\mathrm{5}\times {\mathrm{10}}^{\mathrm{5}}$ m^{−1} to form the initial profile. At the same time, the total scalar concentration was maintained by two Neumann boundary conditions with surface emission equal to 0.0001 and a sink at the domain top equal to $\mathrm{0.0001}\times (\mathrm{1}{\mathit{\lambda}}_{p,\mathrm{0}})$, where λ_{p,0} indicates urban density at the surface.
The detailed design of building arrays can be found in Lu et al. (2023a) and is briefly explained here: while maintaining the same volume of roughness elements (i.e., same mean building height, H_{mean}=16 m), we consider three height distributions, i.e., ${H}_{\mathrm{std}}\in [\mathrm{0}\phantom{\rule{0.125em}{0ex}}\mathrm{m},\phantom{\rule{0.125em}{0ex}}\mathrm{2.8}\phantom{\rule{0.125em}{0ex}}\mathrm{m},\phantom{\rule{0.125em}{0ex}}\mathrm{5.6}\phantom{\rule{0.125em}{0ex}}\mathrm{m}]$, which result in different maximum building height H_{max}. Additionally, we consider three types of urban arrays: (a) continuous distribution of building height with a relatively low H_{max}; (b) clustered configuration with two clusters (three buildings) of low and tall building blocks; and (c) highrise configuration, featuring one tall and a cluster (three) of low building blocks. Each configuration is then run with seven urban packing densities, i.e., ${\mathit{\lambda}}_{p}\in [\mathrm{0.0625},\mathrm{0.1111},\mathrm{0.16},\mathrm{0.25},\mathrm{0.35},\mathrm{0.4444},\mathrm{0.64}]$, which yields a total of 49 simulation cases. Details of the configurations considered are shown in Table 1.
In mesoscale climate models, flow through canopies is often simulated at scales much larger than the typical surface processes, such as turbulent eddies and obstacle wakes. A common approach to derive parameterizations is to apply a time averaging over an interval longer than the scale of the slowest eddies and a spatial averaging over a length larger than the typical spatial deviations in the flow (Raupach and Shaw, 1982; Finnigan, 1985). For momentum, the Reynolds decomposition is first applied to the threedimensional instantaneous equations that decompose mean flow quantities ($\stackrel{\mathrm{\u203e}}{\mathit{\varphi}}$) from their fluctuating components (ϕ^{′}, time or ensemble averaging, $\mathit{\varphi}=\stackrel{\mathrm{\u203e}}{\mathit{\varphi}}+{\mathit{\varphi}}^{\prime}$). Then, flow properties are spatially averaged to match a grid cell of a mesoscale model (horizontal averaging, $\stackrel{\mathrm{\u203e}}{\mathit{\varphi}}=\langle \stackrel{\mathrm{\u203e}}{\mathit{\varphi}}\rangle +\stackrel{\mathrm{\u0303}}{\mathit{\varphi}}$, where $\stackrel{\mathrm{\u0303}}{\mathit{\varphi}}$ is termed a dispersive component).
In urban canopies where the volume fraction occupied by obstacles is not negligible, the intrinsic average (Schmid et al., 2019) that excludes the roughness fraction in the spatial averaging is more favorable given that additional terms compensating for the volume fraction occupied by roughness elements in the canopy are included Blunn et al. (2022). The friction Reynolds number ($Re={u}_{\mathit{\tau}}{H}_{T}/\mathit{\nu}$, where u_{τ}∼0.21 m s^{−1} is the friction velocity, H_{T} is the domain height, and ν is the kinematic viscosity) of the flow dataset is large enough (10^{6}) to neglect viscous stress on momentum and molecular flux on scalars. The resulting onedimensional representation of the momentum and passive scalar equation assuming horizontal homogeneity, negligible Coriolis effect, and buoyancy force reads
where ρ is the air density, u and w are the streamwise and vertical velocity so that $\langle \stackrel{\mathrm{\u203e}}{{u}^{\prime}{w}^{\prime}}\rangle $ and $\langle \stackrel{\mathrm{\u0303}}{u}\stackrel{\mathrm{\u0303}}{w}\rangle $ form the spatially averaged turbulent and dispersive Reynolds stress, respectively. The first term on the RHS of both equations represents vertical transport events, while the second term is a term that compensates for volume change in the intrinsic averaging, where $\mathit{\gamma}\left(z\right)=\mathrm{1}{\mathit{\lambda}}_{p}\left(z\right)$ is the fluid fraction at the vertical index z. The third term, f_{i} and c_{i}, accounts for external sources or sinks, such as pressure gradient driving the flow and surface emission of pollutants. The fourth term of Eq. (1) represents the pressure drop between windward and leeward building facets (Santiago and Martilli, 2010b), which accounts for momentum sink due to buildinginduced drag rising from spatial averaging over the momentum equation. The multilayer model is designed to relate terms in Eq. (1) to the variation of urban surface, where the vertical diffusion of momentum (Sects. 3.1 and 3.2) and drag (Sect. 3.3) are parameterized.
3.1 Characterization of the eddy diffusivity for different flow statistics
The Reynolds averaging over the momentum equation yields nonlinear terms that must be parameterized to close the equation. The most common approach in parameterizing the Reynolds stress in UCMs is based on K theory (as also seen in models compared in Best and Grimmond, 2015, Lipson et al., 2023) to represent the efficiency of turbulent transports by employing the eddydiffusivity coefficient together with the vertical gradient of the flow properties (for momentum, $\langle {K}_{\mathrm{m}}\rangle =\left(\u2329\stackrel{\mathrm{\u203e}}{{u}^{\prime}{w}^{\prime}}\u232a\right./\u2329\frac{\partial \stackrel{\mathrm{\u203e}}{u}}{\partial z}\u232a$). In a variety of ways specifying the eddy diffusivity, the 1.5order closure (Martilli et al., 2002) assuming the transport strength related to the turbulence intensity is commonly employed for urban canopy flow as follows:
where C_{k} is a model constant, l_{k} is a turbulent length scale, and k is TKE that has to be prognostically obtained during the simulation following Eq. (4).
Similar to the momentum equation (Eq. (1)), the first two terms on the RHS represent the transport of TKE, followed by two terms for shear and dispersive TKE production, respectively. D_{k} represents the source of TKE generated through the interaction with the buildings and the airflow, while the last term $\langle \overline{\mathit{\epsilon}}\rangle $ is the viscous dissipation of TKE that can be parameterized as
where C_{ε} and l_{ε} are the model constant and the dissipation length scale, respectively, and have a relationship with turbulent length scale based on the turbulent viscosity C_{μ} (Santiago and Martilli, 2010b) that measures the relative strength of viscous dissipation over turbulent transport:
To date, most models assume the same value for the diffusion coefficient in the momentum, TKE, and scalar equations, possibly due to the lack of data on TKE and scalar budget terms (note that the explicit evaluation of TKE transport and diffusivity, in particular, is only possible from numerical simulations that resolve thirdorder moments, e.g., LES and direct numerical simulations (DNS)). However, previous LES analyses of urban flow indicated a different transport efficiency of varying flow properties (Nazarian et al., 2020). Accordingly, we distinguish between three diffusion coefficients as
The turbulent fluxes, vertical gradient, and eddy diffusivity of momentum, TKE, and scalars are shown in Fig. 1 for staggereduniform configurations for cubic buildings. The roughness flow over other building shapes such as cuboids (Blunn et al., 2022; Li and BouZeid, 2019), cylinders (Agbaglah and Mavriplis, 2019), and spheres (Ye et al., 2016) exhibits different transport behavior. However, most buildings are a combination of cuboids (Oke et al., 2017), and nonurban roughness flow is not suitable to inform urban parameterization. Therefore, the transport behavior over noncuboid roughness is outside of the current scope. Note that the profiles over other designs with building height variability (not shown here) share similar characteristics. Therefore, we chose the uniformheight layout to better reveal the difference among flow properties. Within the canopy, the turbulent fluxes of momentum and TKE are predominantly negative, reflecting downgradient transports driven by their differing production and destruction mechanisms. The turbulent fluxes of passive scalars are positive due to the surface emission being the only source and exhibiting a minor variation across different densities. The eddy diffusivity of scalars (K_{c}) and momentum (K_{m}) maintains a similar shape as a result of the relatively simpler mechanism in their production, destruction, and transport, which justifies the simplification (K_{c}∼K_{m}) that is consistent with Li and BouZeid (2019). However, the eddy diffusivity of TKE (K_{k}) demonstrates a greater complexity with some negative values at the ground and canopy interface that marks the countergradient transport of TKE indicating two vital local shear production regions: (a) the canopytop shear zone is relatively wellstudied and is induced by clearings of roughness elements (Giometto et al., 2016) and the (b) nearground shear zone that is unique in mediumdense layouts (λ_{p}≤0.25) where flow reversal appears and counteracts with the pressure gradient, inducing a significant shear source of TKE (Lu et al., 2023a). This lack of consideration of the nearground flow features partially explains the underestimation of local TKE seen in Nazarian et al. (2020) and can be mitigated by imposing a heightdependent model constant (e.g., Glazunov et al., 2021), but is beyond the scope of this study that aims to refine the turbulent transport characterization within the whole canopy. Apart from the nearground and canopytop inversion, K_{k} shares a similar shape with K_{m} but presents a significantly larger magnitude, which indicates that models using unified K underestimate the transport of TKE and require closer characterization. The distinct behavior of K_{m} in the densest case (purple lines in Fig. 1c) comes from the positive gradient of streamwise velocity in the middle of the canopy and positive Reynolds stress in the lower canopy that is hard to characterize. Such a high density is not rare, especially in some European and South American cities. However, realistic neighborhoods with extremely high density typically exhibit varied horizontal arrangements, leading to different flow characteristics from idealized building blocks (Lu et al., 2023b).
Following the parameterization of length scales in Nazarian et al. (2020), we average the profile of eddy diffusivities shown in Fig. 1 up to H_{mean} to reveal the integral behavior of transport behaviors across configurations. Figure 2 presents the scatter of the vertically averaged eddy diffusivity for momentum K_{m,T} and TKE K_{k,T} without considering dispersive flux, and these include the dispersive contribution (${K}_{\mathrm{m},T+D}$ and ${K}_{k,T+D}$). Characterization of dispersive fluxes within the Ktheory framework enhances the exchange rate of flow properties, which agrees with previous studies such as Akinlabi et al. (2022) where the dispersive flux is mostly downward.
The magnitude of eddy diffusivity for TKE is ∼3.5 times higher than that for momentum (note the different vertical ranges for momentum and TKE in Fig. 2). This considerable difference demonstrates a significant contrast to the parameterization strategy in Nazarian et al. (2020) that assumed ${K}_{\mathrm{m},T+D}={K}_{k,T+D}$, which results in underestimations of TKE within the canopy. In addition, unlike ${K}_{\mathrm{m},T+D}$, which generally decreases with urban density before the layout becomes very dense, ${K}_{k,T+D}$ exhibits a less decisive trend to density. The magnitude difference and variations over densities between ${K}_{\mathrm{m},T+D}$, and ${K}_{k,T+D}$ urge for a distinction in their parameterizations.
The eddy diffusivity was not directly parameterized under the 1.5order closure framework developed in Martilli et al. (2002) but as a turbulent length scale correlated with the dissipation length scale in Eq. (6) (i.e., ${C}_{k}{l}_{k}={C}_{\mathit{\mu}}\frac{{l}_{\mathit{\epsilon}}}{{C}_{\mathit{\epsilon}}}$). Here we follow the same approach as in Nazarian et al. (2020) to parameterize the turbulent viscosity (C_{μ}) and dissipation length scale ($\frac{{l}_{\mathit{\epsilon}}}{{C}_{\mathit{\epsilon}}}$) to obtain the turbulent length scale (C_{k}l_{k}). Considering that the parameterization for dissipation length scale is the same across the scenarios considered in Table 2, the difference in the transport efficiency demonstrated in Sects. 3.1 and 3.2 can be expressed as the difference in the turbulent viscosity (C_{μ}), as shown in Fig. 3.
Figure 3 shows the turbulent viscosity calculated in 49 idealized urban arrays binned over densities for better presentation, where the original parameterization with seven uniformheight cases (Nazarian et al., 2020) is shown in green lines. Generally, C_{μ} increases with density over sparse and medium layouts, which indicates the overweighting trend of viscous dissipation of TKE over the turbulent transport that is consistent with the transition of the flow regime (Oke et al., 2017). Meanwhile, denser layouts amplify the impact of building height variability, leading to a higher standard deviation among binned scenarios, which is consistent with Lu et al. (2023a). Turbulent viscosity for momentum lumping the dispersive flux ${C}_{\mathit{\mu},{K}_{\mathrm{m},T+D}}$ shows a good agreement with that from Nazarian et al. (2020) except for the peaking at higher density. Due to the separate explicit characterization of the buildinginduced dispersive momentum transport, ${C}_{\mathit{\mu},{K}_{\mathrm{m},T+S}}$ is slightly lower than ${C}_{\mathit{\mu},{K}_{\mathrm{m},T+D}}$ but shows a similar shape. The turbulent viscosity of TKE (${C}_{\mathit{\mu},{K}_{k,T+D}}$) is about 3.5 times higher than that of momentum (note the difference on the yaxis limit) and shows exponential growth over density.
3.2 Characterization of the dispersive momentum transport
The momentum transport in urban canopy flow is also modulated by its spatial variability (Eq. 1) induced by the rigid volume of roughness elements that is hard to characterize due to its extreme heterogeneity (Poggi and Katul, 2008b). Nevertheless, flow heterogeneity is not unique to the urban boundary layer but is prevalent in the realworld atmosphere layer due to the surface heating and has been well characterized in PBL parameterizations (e.g., Han and Bretherton, 2019). In the context of convective boundary layer parameterization, the eddy diffusivity (ED) approach based on K theory (as discussed in Sect. 3.1) has been successful in representing neutral and surface boundary layers that promote more homogeneous and isotropic turbulence (Fig. 4a). On the other hand, the flow heterogeneity induced by convection from thermal forcing and moisture phase change has been successfully captured by the massflux (MF) approach (Siebesma et al., 2007). In the effort to model both atmospheric conditions coexisting in PBL, an eddydiffusivity massflux (EDMF) scheme has been developed by partitioning flow fluctuations into contributions from local turbulent mixing (ED) and nonlocal coherent rising and descending air parcels (MF, Lopez‐Gomez et al., 2020; Lu, 2019).
This section explores this alternative approach to characterize the dispersive momentum flux, other than the pragmatic approach developed in Nazarian et al. (2020); Simónmoral et al. (2016), by considering dispersive transports as massflux components under the framework of the EDMF scheme. The dispersive transports of momentum (Fig. 4b) in the urban boundary layer are induced by coherent motions near the building facets and present a strong local distribution that resembles the patchy coherent structures induced by thermal forcing in PBL. Therefore, the similarity in flow heterogeneity between UCL and PBL implies a phenomenological analogy, which favors the introduction of the EDMF approach to multilayer UCMs. Following the formulation in Siebesma et al. (2007), the dispersive flux of momentum is firstly decomposed into individual plumes (a_{i}) enclosing coherent motions with strength $\left(\stackrel{\mathrm{\u203e}}{{u}_{i}}\langle \stackrel{\mathrm{\u203e}}{u}\rangle \right)\left(\stackrel{\mathrm{\u203e}}{{w}_{i}}\langle \stackrel{\mathrm{\u203e}}{w}\rangle \right)$ and then sorted into upward (positive, with a fraction of a_{u}) and downward (negative, with a fraction of a_{d}) contributions. The dispersive flux of momentum then reads
Note that, for simplicity, the entrainment process (Abdella and Petersen, 2000) that represents scalar transport between subdomains (between dispersive structures and the turbulent environment) is not considered here. Using the massflux approach to characterize dispersive momentum flux requires the fraction (${a}_{\mathrm{u}}+{a}_{\mathrm{d}}=\mathrm{1}$) and strength ($\langle \stackrel{\mathrm{\u0303}}{u}\stackrel{\mathrm{\u0303}}{w}{\rangle}_{\mathrm{u}}$, $\langle \stackrel{\mathrm{\u0303}}{u}\stackrel{\mathrm{\u0303}}{w}{\rangle}_{\mathrm{d}}$) of those coherent flow motions to be parameterized based on the variation of the urban arrangement. An example of the spatial distribution of turbulent and dispersive momentum flux is shown in Fig. 4 for the ShighriseSD56 configuration of density λ_{p}=0.25 sampled at z=15 m.
Unlike the turbulent momentum flux, which shows a uniform downward transport with an indecisive correlation with the local roughness elements, a portion of the dispersive counterpart is directly connected to the local roughness, whereas the rest remains relatively homogeneous. Therefore, a proper sampling filter should be applied to enclose these regions before gathering flow properties for parameterization. Again, we adopt a common sampling approach in the development of the EDMF (Couvreux et al., 2010; Sušelj et al., 2012) based on the joint cumulative probability density function (CDF) of the dispersive velocities ($\stackrel{\mathrm{\u0303}}{u},\stackrel{\mathrm{\u0303}}{w}$).
The CDF can be adjusted from [0,1] that corresponds to different integration ranges of dispersive streamwise and vertical velocity (Eq. 9). By lowering the CDF, the sampled regions are adjusted smaller only to enclose regions with outstanding strength. We tested a range of cutoff values of $\mathrm{CDF}\in [\mathrm{0.3},\mathrm{0.4},\mathrm{0.5},\mathrm{0.6},\mathrm{0.7}]$ (test results are shown in the Appendix and in the Supplement), and a cutoff value here is determined as 0.5 for all cases considered. Figure 4b shows the sampled upward and downward transports enclosed by red and blue lines, respectively.
The sampling filter aims to capture the buildinginduced component of the dispersive momentum flux that K theory cannot characterize. Here, we test the correlation between dispersive streamwise velocity ($\stackrel{\mathrm{\u0303}}{u}$) and vertical velocity ($\stackrel{\mathrm{\u0303}}{w}$) to show how the urban structures modify the statistical distribution of dispersive fluxes and the corresponding characterization fulfilled by the sampling. Figure 5 shows the 2D CDF distribution for ShighriseSD56 configuration of density λ_{p}=0.25. The distribution exhibits a more eccentric behavior along the vertical range, with

One tail growing stronger and thinner at the third quadrant ($\stackrel{\mathrm{\u0303}}{u}<\mathrm{0},\stackrel{\mathrm{\u0303}}{w}<\mathrm{0}$) representing the blockage effects of the frontal area of the roughness (Fig. 4d), which is wellcaptured by the sampling strategy (red squares in Fig. 5). The height distribution of building blocks can explain the thinner trend and reflect how the dispersive flux responds to the urban vertical structure.

A relatively symmetric elongated region at the second ($\stackrel{\mathrm{\u0303}}{u}<\mathrm{0},\stackrel{\mathrm{\u0303}}{w}>\mathrm{0}$) and fourth quadrant ($\stackrel{\mathrm{\u0303}}{u}>\mathrm{0},\stackrel{\mathrm{\u0303}}{w}<\mathrm{0}$) representing the leeward flow cavity (Fig. 4b, characterized by the sampled flux at most heights) and flow divergence at the lateral facets (represented by the sampled flux at 2.5 m), respectively. The elongation of both sides reflects a more irregular distribution of downward dispersive flux with heights, which is consistent with the analysis in Lu et al. (2023a).
After subtracting dispersive structures from the sampled regions, the rest of the relatively homogeneous regions may safely resort back to the lumping approach and update Eq. (8) with Eq. (10):
We speculate the reason that dispersive flux has long been ill represented in UCMs is attributed to the inappropriate scaling with the turbulent flux (e.g., Poggi and Katul, 2008b), which usually yields a small fraction without recognizing these two processes are controlled by very different physical processes. Alternatively, the dispersive perturbations are caused by the pressure fields associated with the topography or the changed resistance to flow through the canopy (Finnigan et al., 2015; Finnigan and Shaw, 2008), which implies a possibility of scaling its strength to the external pressure gradient. This was only made possible recently due to the increased computational power and awareness of intrinsic urban heterogeneity (Boeing, 2019; Lu et al., 2023b). Therefore, we hypothesize the strength of sampled flow motions can be parameterized by pressure gradient and urban morphological parameters as Eq. (11):
To gain a direct impression of how the sampled dispersive motions scale with the pressure gradient, we present the profile of the vertical gradient of the sum of two sampled terms from Eq. (10) in Fig. 6. The variation of the dispersive momentum transport exhibits a strong dependence on the flow regime but decouples with the turbulence level (not shown here; see Figs. 9–11 from Lu et al. (2023b) for spatial pattern of TKE as well as turbulent and dispersive momentum flux). The sampled dispersive motion is of a similar magnitude to the only driver in our experiments – pressure gradient collectively leaving such a term unattended is equivalent to only applying half of the pressure gradient to the flow and inevitably degrades model predictability. Considering that realistic urban geometries will exhibit increased complexity and dispersive flux (Giometto et al., 2016), explicit parameterization is necessary. The little variation in the profiles for each density in Fig. 6 indicates building height distribution does not induce significant variance over the sampled dispersive flux among cases within the H_{mean}. Instead, urban density dictates the variation of the sampled flux, whose profile is different between sparse (the first row in Fig. 6) and denser (the second row in Fig. 6) layouts: for the three sparse layouts, a relatively linear increasing trend crossing zero with heights is developed. Denser layouts lead to strong overlapping of building wakes, which helps to maintain a constant strength profile of similar magnitude to the pressure gradient, which indicates a trend of a skimming flow regime (Oke et al., 2017). Above the H_{mean}, sampled fluxes quickly converge to zero except for the highrise configurations over denser layouts (λ_{p}≤0.25) and at the H_{max} of each configuration, demonstrating a minor magnitude but complex pattern across configurations. Therefore, for simplicity, we only provide a parameterization for the sampled fluxes below H_{mean} that scales with the external pressure gradient $\frac{\mathrm{d}p}{\mathrm{d}x}$ as follows:
The above parameterization retains the physical significance of upward and downward dispersive motion demonstrated in Fig. 6. By definition, the scaling coefficients to the external pressure gradient are designed to be lowered over denser layouts, consistent with previous findings where dispersive momentum flux is more significant over sparse arrangements (e.g., Lu et al., 2023b).
3.3 Drag and dissipation parameterization
The presence of roughness elements on the urban surface implicitly modifies turbulent length scales (Li et al., 2020) that regulate the transport mechanisms of flow properties. Explicitly, the buildings block and divert the flow, inducing pressure deficits between their windward and leeward facets that collectively serve as a sink for momentum (Eq. 1) and a source for TKE (Eq. 4). In the following, we revisit the modeling strategy for drag and dissipation in Nazarian et al. (2020) and add the necessary adjustments to reflect the building height variability. By assuming the strength of building effects is proportional to the area facing the wind per cubic meter of outdoor air volume S(z) and the spatially averaged mean wind speed $\langle \stackrel{\mathrm{\u203e}}{U}\rangle $, and a drag coefficient C_{d} encodes the variation of building arrangements, the drag parameterization reads
Accordingly, characterizing drag impacts converts to evaluating the drag coefficient C_{d} that comes in different forms in the literature (e.g., Santiago et al., 2013). The equivalent drag coefficient C_{deq} evaluated following Santiago and Martilli (2010b) requires fewer inputs but maintains a similar performance, as demonstrated by its capacity in MLUCM v2.0 (Nazarian et al., 2020). It will be retained in the present study. Using C_{deq} for the evaluation assumes a constant profile of drag that ensures the integral of the drag over the whole urban canopy is equivalent to that evaluated from LES simulations:
where $\mathrm{\Delta}\langle \stackrel{\mathrm{\u203e}}{p\left(z\right)}\rangle ={p}_{\mathrm{f}}\left(z\right){p}_{\mathrm{b}}\left(z\right)$ is the pressure deficit between the windward and leeward facets of each building block. Due to the low variance of C_{deq} over the impact of height variability except for denser layouts λ_{p}≥0.35, we binned all seven cases over each density in Table 2 to show their mean and standard deviation in Fig. 7 together with that from Nazarian et al. (2020) for uniformstaggered arrays. The excellent agreement with C_{deq} in the present study indicates the parameterization in Nazarian et al. (2020) presents great resilience to the height variability, especially over the sparse layouts. Therefore, updating the constants in the original parameterization requires only a minor change. However, height variability not only extends the vertical range of the urban canopy but also complicates the vertical structure of the drag coefficient and requires extra characterization. Figure 7b shows the ratio between the drag coefficient above the mean building height H_{mean} and C_{deq} that covers a range that is generally larger than unity and shows a decreasing trend with densities, which can be captured by considering the morphological condition of the urban surface from Lu et al. (2023a).
The last term parameterized is the dissipation length scale from Eq. (5). As indicated in the left panel of Fig. 8, the impact of height variability ${l}_{\mathit{\epsilon}}/{C}_{\mathit{\epsilon}}$ below the mean building height only demonstrates a very minor variation for each density range. The variation is more evident above H_{mean}, with a clear shape following the sectional density, and was discussed in Lu et al. (2023a). Despite the fact that applying a sophisticated dissipation and drag parameterization might lead to better model performance (Glazunov et al., 2021; Castro, 2017; Sützl et al., 2021a), it is beyond the scope of the present study, which is focused on refining the transport characterization of UCM. Therefore, we decided to keep the original parameterization from MLUCM v2.0 (Nazarian et al., 2020). For completeness, we briefly describe the parameterization strategy: by segmenting the urban canopy flow into the mixing layer, transition region, and wallbounded layer (Coceal and Belcher, 2004), the parameterization reads
where α_{1}=5.5 is a revised value based on 49 simulations in the present study, and $d={H}_{\mathrm{mean}}{\mathit{\lambda}}_{p}^{\mathrm{0.15}}$ is the displacement height parameterized following Krayenhoff et al. (2015). Equation (15a) depicts a mixing layer type of flow with a constant dissipation efficiency and Eq. (15b) represents the transition to a fully wallbounded flow in Eq. (15c) where the dissipation efficiency decreases with height based on d_{2} and α_{2} parameterized as follows:
In this section, we evaluate the predictability of two modifications to the original MLUCM developed in Nazarian et al. (2020), where the combination of changes proposed in Sects. 3.1 and 3.2 and the original model form four scenarios, as shown in Table 2.
Figure 9 shows the vertical profiles of horizontally averaged velocity, TKE, and turbulent momentum flux calculated with the onedimensional models with different transport characterization and compared with the corresponding LES results for the highrise configuration with H_{std}=5.6 m. The profile comparison for other configurations can be found in the Supplement. We select sparse (λ_{p}=0.1111), medium (λ_{p}=0.25), and dense (λ_{p}=0.4444) layouts to reflect different flow regimes for comparison.
The prediction of streamwise velocity does not present a great variation, where 1DK_{k}EDMF exhibits the highest value. This is explained by the newly introduced massflux term as an explicit source of momentum and differentiation of K_{k} from K_{m} enhancing the momentum transport from the atmosphere to the canopy.
For TKE estimation, 1DK_{k} and 1DK_{k}EDMF with differentiation of K_{k} from K_{m} (Table 2) demonstrate a significant improvement on the TKE profile for medium and dense layouts below H_{mean}. The introduction of EDMF is phenomenologically equivalent to imposing a stronger pressure gradient, which leads to higher velocity and explains the better agreement with LES for the 1DEDMF and 1DK_{k}EDMF scenarios. As a result, 1DK_{k}EDMF yields an even higher TKE due to the introduction of a massflux term to parameterize dispersive momentum flux as an explicit source leading to higher TKE production. The correction of TKE not only leads to better modeling of the velocity profile but also enhances the prediction of urban fluxes such as the overestimation of daytime air temperature by correcting the vertical exchange of heat (Krayenhoff, 2014).
The variation in the performance of scenarios is not pronounced across different height variability, which further consolidates the binning strategy for parameterization in the present study. To evaluate the comprehensive performance of scenarios over different height configurations, Fig. 10 presents the root mean square error (RMSE) between vertical profiles (ranging from $z=[\mathrm{0}\mathrm{3}{H}_{\mathrm{mean}}]$) of LES and four MLUCM scenarios averaged over seven configurations (one uniform and six variable height, as shown in Table 1) over seven densities.
The new parameterizations with refined transport characterization represent an overall improvement over the previous multilayer model (Nazarian et al., 2020). All three scenarios with modifications lead to better predictability for velocity except for the two sparse layouts. The predictability for TKE is substantially better from the 1DK_{k} and 1DK_{k}EDMF scenarios, except for two sparse layouts that can be further improved with a better characterization of the TKE budget terms such as dissipation and wake production (Blunn et al., 2022). The turbulent momentum flux is diagnostically obtained in MLUCM based on the k−l closure ($\langle \stackrel{\mathrm{\u203e}}{{u}^{\prime}{w}^{\prime}}\rangle ={C}_{k}{l}_{k}\langle \stackrel{\mathrm{\u203e}}{k}{\rangle}^{\mathrm{1}/\mathrm{2}}\frac{\partial \u2329\stackrel{\mathrm{\u203e}}{u}\u232a}{\partial z}$) that is highly related to streamwise velocity and TKE. However, it exhibits a consistently contrasting performance trend compared to TKE, which leads to less improvement in denser layouts. This pattern can be attributed to both improvements (Table 2) enhancing the vertical exchange of flow properties within the canopy, which results in a smaller vertical gradient of streamwise velocity. Having the most remarkable performance, 1DK_{k}EDMF with both modifications from Sects. 3.1 and 3.2 is desirable for MLUCM v3.0.
This study refined the characterization of the transport of flow properties for the multilayer urban canopy model via a separate parameterization of TKE diffusion and introduced a “massflux” term for the dispersive transport of momentum. The updated multilayer model demonstrates an improved performance by correcting the underestimation of turbulent exchange and velocity below the mean building height H_{mean}.
By analyzing 49 LES simulations over staggered urban arrays of uniform and variable building height and comprehensive density coverage, we found the turbulent exchange rate of momentum is similar to that for scalars but 3.5 times lower than TKE. However, the previous model (MLUCM v2.0, Nazarian et al., 2020) assumed a unified transport efficiency for momentum and TKE, which curbs the transport of the flow above H_{mean} into the canopy and causes an underestimation of incanopy TKE and wind speed. The accurate characterization of TKE transport into the lower canopy becomes critical as it also controls the turbulent exchange of heat and moisture between buildings and vegetation.
Regarding applications, the enhanced model, which addresses the underprediction of incanopy TKE, has the potential to alleviate the overestimation of daytime air temperature and diurnal temperature range (Krayenhoff et al., 2020). Furthermore, the improved prediction of urban fluxes reinforces the advantages of MLUCM over simpler bulk schemes, emphasizing their ability to realistically simulate vertical diffusion processes in mesoscale atmospheric models (Hendricks et al., 2020). The proposed refinements are not exclusively constrained to urban applications but are also applicable to broader turbulence modeling. For example, the explicit parameterization of dispersive momentum flux could benefit turbulence modeling over other types of roughness elements such as coral reefs (Davis et al., 2021) and vegetative canopies (Finnigan, 2000). In parallel, emphasizing the distinction between TKE transport efficiency and its momentum counterpart underscores the importance of gathering higher flow moments for modeling general atmospheric flow (Ghonima et al., 2017).
We also revealed the distinct horizontal distribution of dispersive momentum flux, which is induced by flow heterogeneity responding to windward (quadrant 3, upward), lateral (quadrant 4, downward), and leeward (quadrant 2, downward) flow patterns. In response to this complication, we applied a sampling filter to the cumulative density function of the dispersive transport to segment the buildinginduced contribution (MF) that scales well with the pressure gradient and the rest of the relatively undisturbed component that can be securely lumped into the turbulent counterpart in the 1.5order turbulent closure (ED). The EDMF framework, having successfully modeled the planetary boundary layer with flow heterogeneity induced by thermal effects, was demonstrated to be also favorable for the representation of flow heterogeneity caused by rigid building volumes. The new framework enables an explicit consideration of the impact of flow heterogeneity in MLUCM, improving our ability to calibrate the model for simulations of complex canopy flow.
With adaptations to the parameterizations of drag coefficient and dissipation length scale to account for building height variability, the model configuration that includes both the updated eddy diffusivity of TKE (K_{k}) and the EDMF scheme yields an improved performance relative to MLUCM v2.0. However, analysis of the flow profiles and length scales indicates future investigation is needed to address the following aspects:

The turbulent length scales show a great variance over heights for layouts with nonuniform building height. A heightdependent parameterization of turbulent length scale and drag profile is necessary to correctly simulate features of urban flow such as the logarithmic or exponential behavior of the wind speed profile (Castro, 2017).

The dispersive transport of TKE and scalars can be considered similarly following the dispersive momentum flux parameterization. In the original implementation of EDMF in PBL parameterizations (Siebesma et al., 2007), the sampled air parcel also exchanges flow properties with the turbulent environment and other air parcels, which could be considered in future developments.

Staggered building arrays were found to be less representative of flow over realistic dense urban neighborhoods (Lu et al., 2023b). The wind direction (Santiago et al., 2013) and thermal effects (Simónmoral et al., 2016) also affect dynamic urban flow parameterizations. Further simulation for model calibration should account for these factors for a more realistic flow characterization.

Similar to the scarcity of evaluations of the eddy diffusivity of flow properties other than momentum, the variation of the TKE budget terms is still poorly considered in urban flow models. Accordingly, a more comprehensive analysis that addresses the relationships between the urban morphology and the TKE budget may provide further improvements.

The applicability of the 1.5order k−l turbulence closure in MLUCM is restricted to isotropic turbulence flow, where the momentum exchange rate is directly correlated with turbulence strength. However, this scaling becomes questionable even in neutral conditions and diminishes under thermal stratification, where turbulence tends to be highly anisotropic (Stiperski and Calaf, 2018). In response to these limitations, employing a higherorder closure that diagnostically assesses momentum transfer (Wilson, 1988; Ayotte et al., 1999) or considering alternative scaling approaches for poorly performing regions (Sun et al., 2020) may contribute to the development of more effective urban canopy models in the future.
The sampling procedure to extract the components of buildinginduced dispersive fluxes was conducted as follows: (a) separating the dispersive fluxes into gradient (downward) and countergradient (upward) components, (b) sorting these two fluxes into a descending order in terms of the magnitude, (c) lowering the cumulative density function of the two fluxes from the largest allows a smaller subdomain to be sampled until subdomains only enclose transport events connected to buildings, and (d) spatially averaging the sampled dispersive fluxes over subdomains for both gradient and countergradient dispersive momentum flux. Figure A1 shows a range of cutoff values ($\mathrm{co}\in [\mathrm{0.3},\mathrm{0.4},\mathrm{0.5},\mathrm{0.6},\mathrm{0.7}]$) for CDF over an example for λ_{p}=0.25. The cutoff value here is determined as 0.5 for all cases considered (an example is provided in Fig. 4b, where blue (red) lines enclose the sampled downward (upward) transports).
The threedimensional distribution of the sampled dispersive structures is shown in Fig. A2 ShighriseSD56 for ${\mathit{\lambda}}_{p}\in [\mathrm{0.0625},\mathrm{0.4444}]$ where the dispersive momentum flux concentrates in the immediate vicinity downstream of the cubes. Structures near the ground contributed significantly to negative stresses (green structures) where the recirculating flows dominate. In response to this distribution (note the dispersive velocity vectors in Fig. 4d), the parameterization for the vertical variation of sampled dispersive fluxes is stronger over the lower part of the canopy over sparse layouts, which is in general agreement with Poggi and Katul (2008a). In the upper levels of the canopy, the countergradient stresses (red structures) are more pronounced at the windward facets due to flow reflection.
The source code and the supporting data of three scenarios for the 1D Multilayer Urban Canopy Model are publicly available at https://doi.org/10.5281/zenodo.10207052 (Lu et al., 2023c) under GPL 3.0 license.
The dataset is available in the Supplement.
The supplement related to this article is available online at: https://doi.org/10.5194/gmd1725252024supplement.
JL, NN, MAH, ESK, and AM collectively developed and planned the study. JL ran the largeeddy simulations designed by NN and modified the model with the help of NN and AM. JL carried out the result analyses and wrote the paper with significant input and critical feedback from NN, MAH, ESK, and AM.
The contact author has declared that none of the authors has any competing interests.
Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.
Simulations were undertaken with the assistance of resources and services from the Australian Government's National Collaborative Research Infrastructure Strategy (NCRIS), with access to computational resources provided by the National Computational Infrastructure (NCI), which is supported by the Australian Government through the National Computational Merit Allocation Scheme. We also thank the anonymous reviewers for their constructive and valuable comments.
This research has been supported by the Australian Research Council Centre of Excellence for Climate Extremes (Grant CE170100023).
This paper was edited by Yongze Song and reviewed by two anonymous referees.
Abdella, K. and Petersen, A.: ThirdOrder Moment Closure Through A MassFlux Approach, Bound.Lay. Meteorol., 95, 303–318, https://doi.org/10.1023/A:1002629010090, 2000. a
Agbaglah, G. and Mavriplis, C.: Threedimensional wakes behind cylinders of square and circular crosssection: early and longtime dynamics, J. Fluid Mech., 870, 419–432, https://doi.org/10.1017/jfm.2019.265, 2019. a
Akinlabi, E., Maronga, B., Giometto, M. G., and Li, D.: Dispersive Fluxes Within and Over a Real Urban Canopy: A LargeEddy Simulation Study, Bound.Lay. Meteorol., 185, 93–128, https://doi.org/10.1007/s10546022007256, 2022. a, b
Ayotte, K. W., Finnigan, J. J., and Raupach, M. R.: A SecondOrder Closure for Neutrally Stratified Vegetative Canopy Flows, Bound.Lay. Meteorol., 90, 189–216, https://doi.org/10.1023/A:1001722609229, 1999. a
Best, M. J. and Grimmond, C. S. B.: Key Conclusions of the First International Urban Land Surface Model Comparison Project, B. Am. Meteorol. Soc., 96, 805–819, https://doi.org/10.1175/BAMSD1400122.1, 2015. a
Blunn, L. P., Coceal, O., Nazarian, N., Barlow, J. F., Plant, R. S., Bohnenstengel, S. I., and Lean, H. W.: Turbulence Characteristics Across a Range of Idealized Urban Canopy Geometries, Bound.Lay. Meteorol., 182, 275–307, https://doi.org/10.1007/s10546021006586, 2022. a, b, c
Boeing, G.: Urban spatial order: street network orientation, configuration, and entropy, Appl. Network Sci., 4, 67, https://doi.org/10.1007/s4110901901891, 2019. a
Bougeault, P. and Lacarrere, P.: Parameterization of OrographyInduced Turbulence in a Mesobeta–Scale Model, Mon. Weather Rev., 117, 1872–1890, https://doi.org/10.1175/15200493(1989)117<1872:POOITI>2.0.CO;2, 1989. a
Brown, M., Lawson, R., DeCroix, D., and Lee, R.: Comparison of Centerline Velocity Measurements Obtained Around 2D and 3D Building Arrays in a Wind Tunnel, article, 1 July 2001, New Mexico, University of North Texas Libraries, UNT Digital Library, https://digital.library.unt.edu/ark:/67531/metadc716934/m1/1/ (last access: 2 April 2024), 2001. a
Castro, I. P.: Are UrbanCanopy Velocity Profiles Exponential?, Bound.Lay. Meteorol., 164, 337–351, https://doi.org/10.1007/s105460170258x, 2017. a, b, c
Cheng, W.C. and PortéAgel, F.: A Simple MixingLength Model for Urban Canopy Flows, Bound.Lay. Meteorol., 181, 1–9, https://doi.org/10.1007/s10546021006500, 2021. a
Christen, A., Rotach, M. W., and Vogt, R.: The Budget of Turbulent Kinetic Energy in the Urban Roughness Sublayer, Bound.Lay. Meteorol., 131, 193–222, https://doi.org/10.1007/s1054600993595, 2009. a
Coceal, O. and Belcher, S. E.: A canopy model of mean winds through urban areas, Q. J. Roy. Meteor. Soc., 130, 1349–1372, https://doi.org/10.1256/qj.03.40, 2004. a
Coceal, O., Thomas, T. G., Castro, I. P., and Belcher, S. E.: Mean Flow and Turbulence Statistics Over Groups of Urbanlike Cubical Obstacles, Bound.Lay. Meteorol., 121, 491–519, https://doi.org/10.1007/s1054600690762, 2006. a
Coceal, O., Dobre, A., Thomas, T. G., and Belcher, S. E.: Structure of turbulent flow over regular arrays of cubical roughness, J. Fluid Mech., 589, 375–409, https://doi.org/10.1017/S002211200700794X, 2007. a
Couvreux, F., Hourdin, F., and Rio, C.: Resolved Versus Parametrized BoundaryLayer Plumes. Part I: A ParametrizationOriented Conditional Sampling in LargeEddy Simulations, Bound.Lay. Meteorol., 134, 441–458, https://doi.org/10.1007/s1054600994565, 2010. a
Davis, K. A., Pawlak, G., and Monismith, S. G.: Turbulence and Coral Reefs, Annu. Rev. Mar. Sci., 13, 343–373, https://doi.org/10.1146/annurevmarine042120071823, 2021. a
Deardorff, J. W.: Stratocumuluscapped mixed layers derived from a threedimensional model, Bound.Lay. Meteorol., 18, 495–527, https://doi.org/10.1007/BF00119502, 1980. a
Finnigan, J.: Turbulence in plant canopies, Annu. Rev. Fluid Mech., 32, 519–571, 2000. a, b
Finnigan, J., Harman, I., Ross, A., and Belcher, S.: Firstorder turbulence closure for modelling complex canopy flows: FirstOrder Turbulence Closure for Canopy Flows, Q. J. Roy. Meteor. Soc., 141, 2907–2916, https://doi.org/10.1002/qj.2577, 2015. a
Finnigan, J. J.: Turbulent Transport in Flexible Plant Canopies, Springer Netherlands, Dordrecht, 443–480, ISBN 9789400953055, https://doi.org/10.1007/9789400953055_28, 1985. a
Finnigan, J. J. and Belcher, S. E.: Flow over a hill covered with a plant canopy, Q. J. Roy. Meteor. Soc., 130, 1–29, https://doi.org/10.1256/qj.02.177, 2004. a
Finnigan, J. J. and Shaw, R. H.: Doubleaveraging methodology and its application to turbulent flow in and above vegetation canopies, Acta Geophys., 56, 534–561, https://doi.org/10.2478/s116000080034x, 2008. a
Frigo, M. and Johnson, S.: FFTW: an adaptive software architecture for the FFT, in: Proceedings of the 1998 IEEE International Conference on Acoustics, Speech and Signal Processing, ICASSP '98 (Cat. No.98CH36181),IEEE, Seattle, WA, USA, ISBN 9780780344280, vol. 3, 1381–1384, https://doi.org/10.1109/ICASSP.1998.681704, 1998. a
Ghonima, M. S., Yang, H., Kim, C. K., Heus, T., and Kleissl, J.: Evaluation of WRF SCM Simulations of Stratocumulus‐Topped Marine and Coastal Boundary Layers and Improvements to Turbulence and Entrainment Parameterizations, J. Adva. Model. Earth Sy., 9, 2635–2653, https://doi.org/10.1002/2017MS001092, 2017. a
Giometto, M. G., Christen, A., Meneveau, C., Fang, J., Krafczyk, M., and Parlange, M. B.: Spatial Characteristics of Roughness Sublayer Mean Flow and Turbulence Over a Realistic Urban Surface, Bound.Lay. Meteorol., 160, 425–452, https://doi.org/10.1007/s1054601601576, 2016. a, b, c
Glazunov, A. V., Debolskiy, A. V., and Mortikov, E. V.: Turbulent Length Scale for Multilayer RANS Model of Urban Canopy and Its Evaluation Based on LargeEddy Simulations, Supercomputing Frontiers and Innovations, 8, https://doi.org/10.14529/jsfi210409, 2021. a, b
Han, J. and Bretherton, C. S.: TKEBased Moist EddyDiffusivity MassFlux (EDMF) Parameterization for Vertical Turbulent Mixing, Weather Forecast., 34, 869–886, https://doi.org/10.1175/WAFD180146.1, 2019. a
Harman, I. N., Böhm, M., Finnigan, J. J., and Hughes, D.: Spatial Variability of the Flow and Turbulence Within a Model Canopy, Bound.Lay. Meteorol., 160, 375–396, https://doi.org/10.1007/s1054601601500, 2016. a
Hendricks, E. A., Knievel, J. C., and Wang, Y.: Addition of Multilayer Urban Canopy Models to a Nonlocal Planetary Boundary Layer Parameterization and Evaluation Using Ideal and Real Cases, J. Appl. Meteorol. Clim., 59, 1369–1392, https://doi.org/10.1175/JAMCD190142.1, 2020. a, b
Kondo, H., Genchi, Y., Kikegawa, Y., Ohashi, Y., Yoshikado, H., and Komiyama, H.: Development of a MultiLayer Urban Canopy Model for the Analysis of Energy Consumption in a Big City: Structure of the Urban Canopy Model and its Basic Performance, Bound.Lay. Meteorol., 116, 395–421, https://doi.org/10.1007/s1054600509055, 2005. a
Krayenhoff, E. S.: A multilayer urban canopy model for neighbourhoods with trees., Ph.D. thesis, University of British Columbia, https://doi.org/10.14288/1.0167084, 2014. a
Krayenhoff, E. S., Santiago, J. L., Martilli, A., Christen, A., and Oke, T. R.: Parametrization of Drag and Turbulence for Urban Neighbourhoods with Trees, Bound.Lay. Meteorol., 156, 157–189, https://doi.org/10.1007/s1054601500286, 2015. a
Krayenhoff, E. S., Jiang, T., Christen, A., Martilli, A., Oke, T. R., Bailey, B. N., Nazarian, N., Voogt, J. A., Giometto, M. G., Stastny, A., and Crawford, B. R.: Urban Climate A multilayer urban canopy meteorological model with trees (BEP Tree): Street tree impacts on pedestrianlevel climate, Urban Clim., 32, 100590, https://doi.org/10.1016/j.uclim.2020.100590, 2020. a, b, c
Kusaka, H., Kondo, H., and Kikegawa, Y.: A simple singlelayer urban canopy model for atmospheric models: comparison with multilayer and slab models, Bound.Lay. Meteorol. 101, 329–358, https://doi.org/10.1023/A:1019207923078, 2001. a
Li, Q. and BouZeid, E.: Contrasts between momentum and scalar transport over very rough surfaces, J. Fluid Mech., 880, 32–58, https://doi.org/10.1017/jfm.2019.687, 2019. a, b, c, d
Li, Q., Bou‐Zeid, E., Grimmond, S., Zilitinkevich, S., and Katul, G.: Revisiting the relation between momentum and scalar roughness lengths of urban surfaces, Q. J. Roy. Meteor. Soc., 146, 3144–3164, https://doi.org/10.1002/qj.3839, 2020. a
Lim, H. D., Hertwig, D., Grylls, T., Gough, H., Reeuwijk, M. v., Grimmond, S., and Vanderwel, C.: Pollutant dispersion by tall buildings: laboratory experiments and LargeEddy Simulation, Exp. Fluids, 63, 92, https://doi.org/10.1007/s00348022034390, 2022. a
Lipson, M. J., Grimmond, S., Best, M., Abramowitz, G., Coutts, A., Tapper, N., Baik, J., Beyers, M., Blunn, L., Boussetta, S., Bou‐Zeid, E., De Kauwe, M. G., De Munck, C., Demuzere, M., Fatichi, S., Fortuniak, K., Han, B., Hendry, M. A., Kikegawa, Y., Kondo, H., Lee, D., Lee, S., Lemonsu, A., Machado, T., Manoli, G., Martilli, A., Masson, V., McNorton, J., Meili, N., Meyer, D., Nice, K. A., Oleson, K. W., Park, S., Roth, M., Schoetter, R., Simón‐Moral, A., Steeneveld, G., Sun, T., Takane, Y., Thatcher, M., Tsiringakis, A., Varentsov, M., Wang, C., Wang, Z., and Pitman, A. J.: Evaluation of 30 urban land surface models in the Urban‐PLUMBER project: Phase 1 results, Q. J. Roy. Meteor. Soc., 150, 126–169, https://doi.org/10.1002/qj.4589, 2023. a, b
Lopez‐Gomez, I., Cohen, Y., He, J., Jaruga, A., and Schneider, T.: A Generalized Mixing Length Closure for Eddy‐Diffusivity Mass‐Flux Schemes of Turbulence and Convection, J. Adv. Model. Earth Sy., 12, e2020MS002161, https://doi.org/10.1029/2020MS002161, 2020. a
Lu, J.: Coherent Structures Sampling and a Stochastic EddyDiffusivity/MassFlux Parameterization for the Transition from Stratocumulus to Cumulus, Master's thesis, UC San Diego, https://escholarship.org/uc/item/6tc324jj (last access: 2 April 2024), 2019. a
Lu, J., Nazarian, N., Anne Hart, M., Krayenhoff, E. S., and Martilli, A.: Representing the effects of building height variability on urban canopy flow, Q. J. Roy. Meteor. Soc., 150, 46–67, https://doi.org/10.1002/qj.4584, 2023a. a, b, c, d, e, f, g, h
Lu, J., Nazarian, N., Hart, M. A., Krayenhoff, E. S., and Martilli, A.: Novel Geometric Parameters for Assessing Flow Over Realistic Versus Idealized Urban Arrays, J. Adv. Model. Earth Sy., 15, e2022MS003287, https://doi.org/10.1029/2022MS003287, 2023b. a, b, c, d, e, f, g, h
Lu, J., Nazarian, N., Anne Hart, M., Krayenhoff, E. S., and Martilli, A: Code for “A onedimensional urban flow model with an Eddydiffusivity Massflux (EDMF) scheme and refined turbulent transport (MLUCM v3.0)” Submitted to GMD, Zenodo [code], https://doi.org/10.5281/zenodo.10207052, 2023c. a
Macdonald, R. W.: Modelling the mean velocity profile in the urban, Bound.Lay. Meteorol. 97, 25–45, https://doi.org/10.1023/A:1002785830512, 2000. a
Maronga, B., Banzhaf, S., Burmeister, C., Esch, T., Forkel, R., Fröhlich, D., Fuka, V., Gehrke, K. F., Geletič, J., Giersch, S., Gronemeier, T., Groß, G., Heldens, W., Hellsten, A., Hoffmann, F., Inagaki, A., Kadasch, E., KananiSühring, F., Ketelsen, K., Khan, B. A., Knigge, C., Knoop, H., Krč, P., Kurppa, M., Maamari, H., Matzarakis, A., Mauder, M., Pallasch, M., Pavlik, D., Pfafferott, J., Resler, J., Rissmann, S., Russo, E., Salim, M., Schrempf, M., Schwenkel, J., Seckmeyer, G., Schubert, S., Sühring, M., von Tils, R., Vollmer, L., Ward, S., Witha, B., Wurps, H., Zeidler, J., and Raasch, S.: Overview of the PALM model system 6.0, Geosci. Model Dev., 13, 1335–1372, https://doi.org/10.5194/gmd1313352020, 2020. a
Martilli, A., Clappier, A., and Rotach, M. W.: An Urban Surface Exchange Parameterisation for Mesoscale Models, Bound.Lay. Meteorol., 104, 261–304, https://doi.org/10.1023/A:1016099921195, 2002. a, b, c
Martilli, A., Santiago, J. L., and Salamanca, F.: On the representation of urban heterogeneities in mesoscale models, Environ. Fluid Mech., 15, 305–328, https://doi.org/10.1007/s1065201393214, 2015. a
McNorton, J. R., Arduini, G., Bousserez, N., Agustí‐Panareda, A., Balsamo, G., Boussetta, S., Choulga, M., Hadade, I., and Hogan, R. J.: An Urban Scheme for the ECMWF Integrated Forecasting System: Single‐Column and Global Offline Application, J. Adv. Model. Earth Sy., 13, e2020MS002375, https://doi.org/10.1029/2020MS002375, 2021. a
Nazarian, N. and Kleissl, J.: Urban Climate CFD simulation of an idealized urban environment: Thermal effects of geometrical characteristics and surface materials, Urban Clim., 12, 141–159, https://doi.org/10.1016/j.uclim.2015.03.002, 2015. a
Nazarian, N., Krayenhoff, E. S., and Martilli, A.: A onedimensional model of turbulent flow through “urban” canopies (MLUCM v2.0): updates based on largeeddy simulation, Geosci. Model Dev., 13, 937–953, https://doi.org/10.5194/gmd139372020, 2020. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x, y, z, aa
Oke, T. R., Mills, G., Christen, A., and Voogt, J. A.: Urban Climates, Cambridge University Press, https://doi.org/10.1017/9781139016476, 2017. a, b, c, d, e
Piacsek, S. A. and Williams, G. P.: Conservation properties of convection difference schemes, J. Comput. Phys., 6, 392–405, https://doi.org/10.1016/00219991(70)900380, 1970. a
Poggi, D. and Katul, G. G.: The effect of canopy roughness density on the constitutive components of the dispersive stresses, Exp. Fluids, 45, 111–121, https://doi.org/10.1007/s0034800804677, 2008a. a
Poggi, D. and Katul, G. G.: Micro and macrodispersive fluxes in canopy flows, Acta Geophys., 56, 778–799, https://doi.org/10.2478/s1160000800297, 2008b. a, b, c
Raupach, M. R. and Shaw, R. H.: Averaging procedures for flow within vegetation canopies, Bound.Lay. Meteorol., 22, 79–90, https://doi.org/10.1007/BF00128057, 1982. a
Redon, E., Lemonsu, A., and Masson, V.: An urban trees parameterization for modeling microclimatic variables and thermal comfort conditions at street level with the Town Energy Balance model (TEBSURFEX v8.0), Geosci. Model Dev., 13, 385–399, https://doi.org/10.5194/gmd133852020, 2020. a
Ross, A. N.: Boundarylayer flow within and above a forest canopy of variable density, Q. J. Roy. Meteor. Soc., 138, 1259–1272, https://doi.org/10.1002/qj.989, 2012. a
Salamanca, F., Krpo, A., and Martilli, A.: A new building energy model coupled with an urban canopy parameterization for urban climate simulations – part I. formulation, verification, and sensitivity analysis of the model, Theor. Appl. Climatol., 99, 331–344, https://doi.org/10.1007/s0070400901429, 2010. a
Santiago, J. L. and Martilli, A.: A Dynamic Urban Canopy Parameterization for Mesoscale Models Based on Computational Fluid Dynamics ReynoldsAveraged NavierStokes Microscale Simulations, Bound.Lay. Meteorol., 137, 417–439, https://doi.org/10.1007/s1054601095384, 2010a. a
Santiago, J. L. and Martilli, A.: A Dynamic Urban Canopy Parameterization for Mesoscale Models Based on Computational Fluid Dynamics ReynoldsAveraged Navier–Stokes Microscale Simulations, Bound.Lay. Meteorol., 137, 417–439, https://doi.org/10.1007/s1054601095384, 2010b. a, b, c, d, e
Santiago, J. L., Coceal, O., and Martilli, A.: How to Parametrize UrbanCanopy Drag to Reproduce WindDirection Effects Within the Canopy, Bound.Lay. Meteorol., 149, 43–63, https://doi.org/10.1007/s105460139833y, 2013. a, b
Santiago, J. L., Krayenhoff, E. S., and Martilli, A.: Urban Climate Flow simulations for simplified urban configurations with microscale distributions of surface thermal forcing, Urban Clim., 9, 115–133, https://doi.org/10.1016/j.uclim.2014.07.008, 2014. a
Schmid, M. F., Lawrence, G. A., Parlange, M. B., and Giometto, M. G.: Volume Averaging for Urban Canopies, Bound.Lay. Meteorol., 173, 349–372, https://doi.org/10.1007/s10546019004703, 2019. a
Schoetter, R., Kwok, Y. T., de Munck, C., Lau, K. K. L., Wong, W. K., and Masson, V.: Multilayer coupling between SURFEXTEBv9.0 and MesoNHv5.3 for modelling the urban climate of highrise cities, Geosci. Model Dev., 13, 5609–5643, https://doi.org/10.5194/gmd1356092020, 2020. a
Siebesma, A. P., Soares, P. M. M., and Teixeira, J.: A Combined EddyDiffusivity MassFlux Approach for the Convective Boundary Layer, J. Atmos. Sci., 64, 1230–1248, https://doi.org/10.1175/JAS3888.1, 2007. a, b, c
Simónmoral, A., Santiago, J. L., and Martilli, A.: Effects of Unstable Thermal Stratification on Vertical Fluxes of Heat and Momentum in Urban Areas, Bound.Lay. Meteorol., 163, 103–121, https://doi.org/10.1007/s1054601602114, 2016. a, b, c, d
Sirko, W., Kashubin, S., Ritter, M., Annkah, A., Bouchareb, Y. S. E., Dauphin, Y., Keysers, D., Neumann, M., Cisse, M., and Quinn, J.: ContinentalScale Building Detection from High Resolution Satellite Imagery, arXiv [preprint], https://doi.org/10.48550/arXiv. 2107.12283, 2021. a
Stiperski, I. and Calaf, M.: Dependence of near‐surface similarity scaling on the anisotropy of atmospheric turbulence, Q. J. Roy. Meteor. Soc., 144, 641–657, https://doi.org/10.1002/qj.3224, 2018. a
Sun, J., Takle, E. S., and Acevedo, O. C.: Understanding Physical Processes Represented by the Monin–Obukhov Bulk Formula for Momentum Transfer, Bound.Lay. Meteorol., 177, 69–95, https://doi.org/10.1007/s10546020005465, 2020. a
Sušelj, K., Teixeira, J., and Matheou, G.: Eddy Diffusivity/Mass Flux and Shallow Cumulus Boundary Layer: An Updraft PDF Multiple Mass Flux Scheme, J. Atmos. Sci., 69, 1513–1533, https://doi.org/10.1175/JASD11090.1, 2012. a
Sützl, B. S., Rooney, G. G., Finnenkoetter, A., Bohnenstengel, S. I., Grimmond, S., and van Reeuwijk, M.: Distributed urban drag parametrization for sub‐kilometre scale numerical weather prediction, Q. J. Roy. Meteor. Soc., 147, 3940–3956, https://doi.org/10.1002/qj.4162, 2021a. a, b
Sützl, B. S., Rooney, G. G., and van Reeuwijk, M.: Drag Distribution in Idealized Heterogeneous Urban Environments, Bound.Lay. Meteorol., 178, 225–248, https://doi.org/10.1007/s10546020005670, 2021b. a
Williamson, J.: Lowstorage RungeKutta schemes, J. Comput. Phys., 35, 48–56, https://doi.org/10.1016/00219991(80)900339, 1980. a
Wilson, J. D.: A secondorder closure model for flow through vegetation, Bound.Lay. Meteorol., 42, 371–392, https://doi.org/10.1007/BF00121591, 1988. a
Xie, Z.T., Coceal, O., and Castro, I. P.: LargeEddy Simulation of Flows over Random Urbanlike Obstacles, Bound.Lay. Meteorol., 129, 1–23, https://doi.org/10.1007/s1054600892901, 2008. a
Ye, Q., Schrijer, F. F., and Scarano, F.: Geometry effect of isolated roughness on boundary layer transition investigated by tomographic PIV, International J. Heat Fluid Flow, 61, 31–44, https://doi.org/10.1016/j.ijheatfluidflow.2016.05.016, 2016. a
Yuan, C., Shan, R., Zhang, Y., Li, X.X., Yin, T., Hang, J., and Norford, L.: Multilayer urban canopy modelling and mapping for traffic pollutant dispersion at high density urban areas, Sci. Total Environ., 647, 255–267, https://doi.org/10.1016/j.scitotenv.2018.07.409, 2019. a
 Abstract
 Introduction
 Methodology
 Assessment of the multilayer urban canopy model (MLUCM)
 MLUCM performance with updated parameterizations
 Summary and conclusions
 Appendix A: Sampling of the dispersive momentum flux
 Code availability
 Data availability
 Author contributions
 Competing interests
 Disclaimer
 Acknowledgements
 Financial support
 Review statement
 References
 Supplement
massfluxterm. These adjustments enhance the model's performance, offering more reliable temperature and surface flux estimates.
 Abstract
 Introduction
 Methodology
 Assessment of the multilayer urban canopy model (MLUCM)
 MLUCM performance with updated parameterizations
 Summary and conclusions
 Appendix A: Sampling of the dispersive momentum flux
 Code availability
 Data availability
 Author contributions
 Competing interests
 Disclaimer
 Acknowledgements
 Financial support
 Review statement
 References
 Supplement