Abstract

The ability to make forecasts depends on atmospheric predictability and the growth of errors. It has recently been shown that the predictability of urban boundary layers differs in important respects from that of the free atmosphere on the mesoscale and larger; in particular, nonlinearity may play a less prominent role in the error evolution. This paper investigates the applicability of linear theory to the error evolution in turbulent street-canyon flow. Using large-eddy simulation, streamwise aspect ratios between 0.15 and 1.50, and identical-twin experiments, it is shown that the growth rate of the error kinetic energy can be estimated from Eulerian averages and that linear theory provides insight into the spatial structure of the error field after saturation. The results should be applicable to cities with deep and closely spaced canyons. Implications for data assimilation and modeling are discussed.

1. Introduction

High-resolution limited-area modeling with horizontal resolution of 1 km or less is now commonplace. Given long-term trends, it seems inevitable that the march toward smaller scales will continue. This progress should offer significantly more detailed predictions within cities, where an increasing percentage of the world’s population lives.

Nevertheless, urban scales pose unique challenges. On the one hand, small-scale urban heterogeneity—notably the presence of buildings, which strongly affect winds and temperatures (e.g., Barlow 2014)—means that building-resolving numerical models are computationally expensive: it is likely that, except for special applications, they will be restricted to modest domains for the foreseeable future. On the other hand, urban observations, especially in the vertical direction, tend to be sparse and the development of a dense observational network, along the lines of the traditional meteorological one, would be a huge undertaking. It is therefore important to consider the potential benefits that could be obtained from high-resolution modeling in general and the assimilation of observations in particular. Such questions fall within the general area of predictability.

There has been surprisingly little work on the predictability of the urban atmosphere, but a recent study with a building-resolving, large-eddy-simulation model of flow over a street canyon (Lo and Ngan 2015) indicates that predictability on the smallest urban scales differs greatly from that of the large-scale atmosphere. Instead of a cascade of errors from small to large scales, that is, a so-called inverse error cascade (e.g., Tribbia and Baumhefner 2004; Ngan et al. 2009), for a street canyon of unit aspect ratio and a neutral atmosphere, the error growth occurs primarily via linear instability. This represents a major departure from homogeneous turbulence, for which the error growth is driven by turbulent cascades, that is, strong nonlinearity (Lorenz 1969; Leith and Kraichnan 1972). In physical terms the “butterfly effect,” which encapsulates the conventional paradigm for atmospheric predictability, is not strictly relevant (cf. Durran and Gingrich 2014). The presence of buildings (or solid boundaries) has profound effects on the predictability of small-scale urban flows.

A previously unexplored implication of Lo and Ngan (2015) is that it may be possible to model the error evolution on urban scales using the linear dynamics. Because standard algorithms such as the extended Kalman filter and incremental four-dimensional variational data assimilation (4D-Var) assume linearity (Daley 1993; Kalnay 2002), this suggests that data assimilation on urban scales may not require new approaches. Research on urban data assimilation has been largely confined to the coupling of mesoscale and urban computational fluid dynamics (CFD) models (Schlunzen et al. 2011; Yamada and Koike 2011; Zajaczkowski et al. 2011) using techniques such as nudging (Stauffer and Seaman 1994). Nevertheless, the development of urban-scale data assimilation would be greatly accelerated if existing systems and approaches, the majority of which assume linearity, could be utilized or leveraged in some way. Despite the considerable interest in nonlinear algorithms such as the ensemble Kalman filter (e.g., Kalnay 2002; Lorenc 2003), linear methods (e.g., Clayton et al. 2013) will probably continue to play a role in research and forecasting.

In this paper the relevance of linear error dynamics to urban flow will be assessed using identical-twin experiments and large-eddy simulation (LES) of turbulent flow in a street canyon. The actual nonlinear growth will be calculated by taking the difference between control and perturbed simulations. The relevance of linear theory after the perturbation has reached statistical equilibrium will be quantified by calculating the spatial correlation between the error field and the (spatially dependent) linear growth rates. A wide range of aspect ratios spanning the isolated-roughness, wake-interference, and skimming-flow regimes (Oke 1988) will be considered because the nature of the error dynamics could change with the flow regime or aspect ratio.

After reviewing the methods for the large-eddy simulations and identical-twin experiments (section 2), results for the error growth (section 3), spatial structure of the errors (section 4), and aspect-ratio dependence (section 5) will be presented. The key finding is that the error dynamics before and after saturation can be usefully modeled with linear theory. The physical significance of the results will be discussed in section 6. In particular, potential applications to data assimilation and modeling will be considered.

2. Methods

a. Large-eddy simulations

Numerical simulations of turbulent flow were performed with the Parallelized LES Model (PALM; Raasch and Schröter 2001; Maronga et al. 2015), which uses the Deardorff subgrid-scale scheme (Deardorff 1980). PALM has been widely used in urban meteorological research (e.g., Letzel et al. 2008, 2012; Inagaki et al. 2012; Park and Baik 2013; Lo and Ngan 2015). Technical details on the numerical schemes may be found in Maronga et al. (2015).

The setup of the simulations closely follows that described in Ngan and Lo (2016). A single street canyon with height H, width W, length L, and streamwise aspect ratio AR ≡ H/W is located at the center of the domain. The street-canyon geometry is very simple, but it can be thought of as the building block of more complicated urban morphologies (cf. Li et al. 2006). The aspect ratio was varied by fixing the separation W and streamwise grid spacing Δx while scaling H and L with AR. To be more precise:

 
formula
 
formula

where Lz and Ly are the vertical and spanwise domain dimensions for AR = 1 and Δz0 and Δy0 are the corresponding grid spacings. The number of grid points in each direction was fixed.

Other prescriptions could be adopted (e.g., the domain size could be kept fixed), but the one described above maintains the number of points within the canyon while keeping Lz/H and ARyL/W fixed. The same approach was used in Ngan and Lo (2016) to verify predictions from the numerical Green’s function.

A standard model configuration was adopted (e.g., Letzel et al. 2008): cyclic boundary conditions in the lateral directions, free-slip boundary conditions at the top, and velocity at solid surfaces parameterized by a roughness length of 0.1 m. The flow was forced by a specified streamwise pressure gradient, and thermal effects were neglected (i.e., the density ρ = 1 kg m−3).

The key parameters are summarized in Table 1. The streamwise grid spacing Δx0 = 0.5 m, implying 100 streamwise grid points within the canyon, while the vertical grid spacing satisfies 0.15 ≤ Δz ≤ 1.5 m. The grid is isotropic only for AR = 0.50, but the anisotropy is fairly weak (Δzx ranges between 0.6 and 2.5). The horizontal resolution exceeds recommendations for urban CFD (Franke et al. 2007), and 20–60 grid points inside the canyon were successfully used in previous street-canyon studies (Liu and Barth 2002; Cui et al. 2004). Neighborhoods composed of several street canyons and horizontal domains of O(2 km × 2 km) have been simulated with coarser grid spacings of approximately 3 m (Hanna et al. 2006; Letzel et al. 2012).

Table 1.

Summary of model parameters for a unit-aspect-ratio canyon; Δx0, Δy0, and Δz0 are the grid spacings, Lx, Ly, and Lz are the domain dimensions, and W is the canyon width.

Summary of model parameters for a unit-aspect-ratio canyon; Δx0, Δy0, and Δz0 are the grid spacings, Lx, Ly, and Lz are the domain dimensions, and W is the canyon width.
Summary of model parameters for a unit-aspect-ratio canyon; Δx0, Δy0, and Δz0 are the grid spacings, Lx, Ly, and Lz are the domain dimensions, and W is the canyon width.

Results for three different aspect ratios (AR = 0.15, 0.60, and 1.50) will be highlighted. These values fall within the isolated-roughness, wake-interference, and skimming-flow regimes for an isolated street canyon (Oke 1988; Hunter et al. 1990). The regime transitions, which occur at approximately AR = 0.3 and AR = 0.7, correspond to changes in the structure of the turbulent wakes; alternatively, they can be associated with changes in long-range vortex interactions between the shear layers at the ground and roof level (Ngan and Lo 2016).1 We shall examine whether the error dynamics depend on the flow regime.

The model has been successfully validated against wind-tunnel measurements of turbulence (Letzel et al. 2008; Lo and Ngan 2015) and pollutant concentration (Lo and Ngan 2017). As the configuration closely follows that of previous studies (Lo and Ngan 2015; Ngan and Lo 2016), the validation is not repeated here.

The predictability results cannot be validated in an engineering sense because appropriate experimental data do not exist. Consequently, the perspective adopted herein is the one that is widely practiced in atmospheric science wherein an imperfect but thoroughly tested model is used to make predictions that cannot be directly validated with experimental data. Nevertheless, the robustness of the results will be investigated.

b. Identical-twin experiments

Following previous studies of the predictability of turbulent flow (e.g., Tribbia and Baumhefner 2004; Ngan et al. 2009; Ngan and Eperon 2012; Lo and Ngan 2015), identical-twin experiments were performed and analyzed. A spinup simulation (of length t0 = 5000 s) is continued with and without small-scale perturbations ε(t0). The perturbed run is referred to as up and the unperturbed or control run is referred to as uc. Errors are diagnosed by differencing the control and perturbation, that is, ε(tt0) ≡ up(tt0) − uc(tt0).

The initial perturbation was applied in spectral space. The perturbation was generated by randomizing the phases of the complex spectral coefficients of uc(t0) and scaling them by a perturbation amplitude α that is nonzero only for streamwise wavenumbers near the Nyquist scale. For the runs reported below, α = 10−4. The results do not depend sensitively on the specification of the perturbation (Lo and Ngan 2015).

The upscale propagation of the error (Leith and Kraichnan 1972) is associated with the inverse error cascade. It leads to a loss of predictability in a deterministic sense, that is, with specific realizations of the perturbation. In principle, loss of deterministic predictability does not preclude statistical predictability of the ensemble, as is widely recognized in long-range forecasting (Palmer and Hagedorn 2006); the distinction between statistical and deterministic predictability is less important for homogeneous turbulence, however.

The preceding formulation does not include model error. The effects of continuous forcing were considered by adding spatially uncorrelated noise to the velocity field at each time step; to avoid numerical problems, the maximum amplitude of the continuous forcing was much smaller (by a factor of 103) than that of the initial error. Continuous forcing will be considered in section 3a.

As in Lo and Ngan (2015), the dimensionless time

 
formula

will be adopted, where the canyon circulation time scale Tc is given by

 
formula

The terms Urms and Wrms denote root-mean-square horizontal and vertical velocities within the canyon, and Tc represents a good estimate of the Lagrangian circulation time of particles: calculations with a Lagrangian stochastic model indicate that the mean residence time of particles is approximately 1.2Tc, with 75% of particles escaping without being reentrained (Lo and Ngan 2017). The nondimensionalization is introduced for convenience and to facilitate interpretation: the results do not depend on the choice of Tc. Growth rates, for example, have been calculated using the dimensional time t.

Note that Tc includes contributions from the mean flow and the turbulence. Inside the canyon, however, the turbulence is fairly weak as the TKE is maximized near the roof level (Liu and Barth 2002; Cui et al. 2004). Consequently, Tc largely characterizes the mean flow. Indeed, TcT1, where T1 is defined analogous to Eq. (4); that is,

 
formula

where U1 and W1 denote spatial averages of the absolute values within the canyon. The difference between T1 and Tc is approximately 10%–15% for AR = 0.15 and 1.50. This suggests that Tc is appropriate for understanding the linear contribution to the error dynamics.

c. Physical-space diagnostics

The errors will be analyzed in physical space. The error kinetic energy or raw error is given by

 
formula

where || || is the Eulerian norm. The normalized error is

 
formula

where E = 0.5||uc||2 is the kinetic energy. There is a complete loss of predictability in a deterministic sense when r is equal to 1.

Three regions will be considered: 1) inside the canyon (z/H < 1); 2) outside the canyon (1 < z/H < Lz/H); and 3) the entire computational domain (z/H < Lz/H). They will be averaged in various ways. The corresponding three-dimensional spatial averages are denoted by a caret; we will focus on averages within the canyon, for example, . Cross sections in the xz plane—for example, 〈ΔE〉(x, z)—are obtained by averaging over y; the angle brackets represent a spanwise average. Time averaging is denoted by the overbar.

The spatial structure of ΔE(x) will be compared with that of other fields. This can be accomplished with spatial correlation coefficients between fields ψ and χ, namely,

 
formula

d. Linear diagnostics

The main concern of this paper is the applicability of linear theory to the error evolution. From the equations of motion, the error obeys (Lo and Ngan 2015)

 
formula

where Ui is the control velocity, D/Dt ≡ ∂/∂t + Ui(∂/∂xi), pc and pp are the control and perturbation pressures, and ν is the (turbulent) viscosity. Note that εi is independent of any externally prescribed forcing. Following the terminology in boundary layer meteorology (e.g., Foken 2008), the terms on the right-hand side represent the contributions of shear production, nonlinearity, pressure transport, and dissipation to the error budget. Vertical profiles of the terms (evaluated at a fixed streamwise location x) are constructed by averaging in the streamwise direction. The error growth is described by the linear equation

 
formula

when the shear production or linear term dominates.

The actual error growth will be compared with that implied by the linear growth rates. For steady, spatially uniform flow, the growth rates Γi are determined by the eigenvalues of −∂Ui/∂xj (Ngan et al. 2004). In more general terms, the growth rates depend on matrix products of the −∂Ui/∂xj evaluated along Lagrangian trajectories. If one ignores the Lagrangian evolution for simplicity, the error in the large-time limit is given by the largest eigenvalue λ:

 
formula

where the tilde represents the space–time average over the canyon; will be a good estimator of the actual canyon-averaged growth rate when nonlinearity is weak and the variation of −∂Ui/∂xj along trajectories is small. Henceforth we shall refer to Eq. (11) as the “linear prediction”; it is also required that the spatial and temporal variability of the velocity gradient be small, however.

In terms of physics, the error will grow for a finite period of time, which is referred to as the saturation time τsat, before reaching a statistically steady value. Time τsat is estimated from the time series of r.2 The growth rate is evaluated from a least squares fit for τ < τsat. The time averaging for the largest eigenvalue is taken over τ > τsat, which represents a more demanding comparison for .

3. Error growth

In this section we examine the perturbation growth for AR = 0.15, 0.60, and 1.50. Table 2 summarizes the key time scales Tc and τsat for these aspect ratios; the error estimates for τsat represent a rough estimate of the interval over which saturation occurs. The aspect-ratio dependence will be analyzed more carefully in section 5.

Table 2.

Summary of statistics for the three aspect ratios studied in sections 3 and 4: Tc is the canyon circulation time scale [Eq. (4)] and τsat is the time required for the perturbation to stop growing (as determined from the change in the local slope of , the canyon-averaged relative error); β is the normalized error difference [Eq. (12)]. The numbers in parentheses are the error ranges for τsat and β.

Summary of statistics for the three aspect ratios studied in sections 3 and 4: Tc is the canyon circulation time scale [Eq. (4)] and τsat is the time required for the perturbation to stop growing (as determined from the change in the local slope of , the canyon-averaged relative error); β is the normalized error difference [Eq. (12)]. The numbers in parentheses are the error ranges for τsat and β.
Summary of statistics for the three aspect ratios studied in sections 3 and 4: Tc is the canyon circulation time scale [Eq. (4)] and τsat is the time required for the perturbation to stop growing (as determined from the change in the local slope of , the canyon-averaged relative error); β is the normalized error difference [Eq. (12)]. The numbers in parentheses are the error ranges for τsat and β.

a. Time series

Figure 1 plots time series of the relative error for the different aspect ratios and averaging regions. There is approximately exponential growth for τ less than approximately 1. Physical growth occurs on the linear time scale. For τ greater than approximately 1, the perturbation saturates; that is, the time average of approaches a constant. The value of is greatest inside the canyon and smallest above it, in agreement with the analysis for AR = 1.0 (Lo and Ngan 2015). The behavior of the raw error is similar. In particular, the scaled growth rates (slopes) of and differ by less than 1%, reflecting the fact that the perturbation has a relatively small effect on the total energy.

Fig. 1.

Time series of the (dimensionless) relative error for AR of (a) 0.15, (b) 0.60, and (c) 1.50. Results for the different averaging regions, i.e., the canyon (green times signs), overlying atmosphere (blue plus signs), and entire computational domain (red circles), are plotted separately. The linear prediction (black triangles) is also shown; the growth rate represents a long-time average over the canyon. It works moderately well inside the canyon.

Fig. 1.

Time series of the (dimensionless) relative error for AR of (a) 0.15, (b) 0.60, and (c) 1.50. Results for the different averaging regions, i.e., the canyon (green times signs), overlying atmosphere (blue plus signs), and entire computational domain (red circles), are plotted separately. The linear prediction (black triangles) is also shown; the growth rate represents a long-time average over the canyon. It works moderately well inside the canyon.

A complete loss of predictability in a deterministic sense occurs when ≥ 1. This happens only for AR = 0.15, because is approximately 1.2 in this case. By contrast, some residual predictability persists for AR = 0.60 and 1.50, because is approximately 0.60 and 0.72, respectively. This agrees with the structure of the flow regimes. In the wake-interference and skimming-flow regimes, large canyon vortices are present and turbulence is inhibited within them.

The agreement between and depends on the aspect ratio. For AR = 0.15, the linear estimate underpredicts the actual growth; for AR = 0.60 and short times τ of less than approximately 0.15, there is excellent agreement; for AR = 1.50, there is fairly good agreement for intermediate times of approximately 0.3 < τ < 0.7. The agreement can be quantified with the normalized difference

 
formula

which ranges between 0.1 and 0.4 (Table 2).

b. Nonlinearity

The departure of the actual error growth from the linear prediction (Fig. 1; Table 2) cannot be attributed to strong nonlinearity. Vertical profiles of the terms in Eq. (9) demonstrate that the nonlinear and pressure terms remain small for τ < τsat. For AR = 0.15 (Fig. 2), the linear term dominates until τ is approximately 0.3; for AR = 1.50 (Fig. 3), the terms become comparable when τ is approximately 1. There are similar results for AR = 0.60 as well as for other vertical profiles (not shown).

Fig. 2.

Vertical profiles of linear (red solid), pressure (blue dashed), and nonlinear (green dot–dashed) error terms [m2 s−3; Eq. (9)] at the downwind wall x = 0.5W for AR = 0.15. In all cases, the linear term dominates prior to saturation of the perturbation; there is similar behavior at the center of the canyon.

Fig. 2.

Vertical profiles of linear (red solid), pressure (blue dashed), and nonlinear (green dot–dashed) error terms [m2 s−3; Eq. (9)] at the downwind wall x = 0.5W for AR = 0.15. In all cases, the linear term dominates prior to saturation of the perturbation; there is similar behavior at the center of the canyon.

Fig. 3.

As in Fig. 2, but at x = 0.5W for AR = 1.50. The behavior is similar to that for AR = 0.15: the linear term dominates prior to saturation.

Fig. 3.

As in Fig. 2, but at x = 0.5W for AR = 1.50. The behavior is similar to that for AR = 0.15: the linear term dominates prior to saturation.

The magnitude of the nonlinear term does not remain small for τ > τsat. Near the roof level, it becomes stronger than the linear term for AR = 0.15 and 1.50. Deep inside the canopy, the terms are roughly comparable.

The relative importance of nonlinearity to linearity is determined by integrating the absolute values of the nonlinear and linear terms from the ground to the roof level. Prior to saturation, the ratio is O(103) for AR = 0.15, O(104) for AR = 0.60, and O(105) for AR = 1.50. Following saturation, the ratio is approximately 0.7–0.8.

c. Lagrangian variability

A likely explanation for the trends described earlier (Fig. 1) is that the assumption of a steady, spatially uniform strain field does not hold equally well for all aspect ratios. While the actual perturbation growth depends on the strain experienced along a trajectory, represents a spatial average. Hence the accuracy of the linear estimate depends on the degree to which the velocity field approximates the special case of a frozen, spatially uniform strain field. Since the structure of the velocity field is determined by the aspect ratio and flow regime, β should also depend on AR. In strict terms, and should be calculated from the cumulative growth along a Lagrangian trajectory, as in the calculation of Lyapunov exponents (e.g., Kalnay 2002). Relatively small differences between Eulerian and Lagrangian averages can be important because the perturbation growth is exponential in .

Averages along Lagrangian trajectories have not been calculated because they require a Lagrangian particle model. Nonetheless, the results of Fig. 1 are consistent with the preceding physical argument. One expects the worst agreement for the isolated-roughness regime (AR = 0.15), for which the turbulence is localized near the walls and the large-scale canyon circulation is relatively weak. One expects improved agreement for the wake-interference (AR = 0.60) and skimming-flow (AR = 1.50) regimes as the flow becomes more spatially coherent and a large-scale circulation develops. Note that β underrepresents the agreement for AR = 1.50 and τ greater than approximately 0.2. For a more spatially coherent velocity field, it may take time for the perturbation to align with the strain axes [cf. Ngan et al. (2004) and references therein].

Representing a turbulent flow with a single parameter (e.g., ) is difficult. In a typical case, analytical predictions for boundary layer turbulence (e.g., the classical log profile or Monin–Obukhov similarity theory) require tuning with semiempirical parameters for good quantitative agreement to be achieved (e.g., Foken 2008). Given that no tuning has been applied, agreement of the actual and predicted growth rates within a factor of 2 seems reasonable.

A similar comparison has not been carried out using large-scale atmospheric winds. One expects, however, that strong nonlinearity prior to saturation should degrade agreement. A closely related calculation for homogeneous, quasi-two-dimensional turbulence (Ngan et al. 2004), in which alignment effects and variability along Lagrangian trajectories play an important role, yields larger errors (β ≈ 0.6).

d. Numerical effects

The behavior described above is not an artifact of resolution or model configuration. Figure 4 considers the sensitivity to horizontal and vertical resolution. Because Δz scales with AR [Eq. (2)], the sensitivity to Δx and Δz also depends on AR. For AR= 0.15 (Fig. 4a), the systematic underprediction of by (Fig. 1a) persists: the agreement worsens as Δx increases toward the value listed in Table 1. For AR = 1.0 (Fig. 4b), the linear estimate works well for small τ and Δz → Δx (Figs. 1b,c). There is analogous behavior when the spanwise domain size Ly is varied (not shown). Figure 5 shows that the results for AR = 0.15 are almost identical for other realizations of the initial perturbation.

Fig. 4.

Resolution tests for (a) AR = 0.15 and Δx and (b) AR = 1.0 and Δz. The linear prediction is shown with symbols. The results are consistent with Fig. 1.

Fig. 4.

Resolution tests for (a) AR = 0.15 and Δx and (b) AR = 1.0 and Δz. The linear prediction is shown with symbols. The results are consistent with Fig. 1.

Fig. 5.

Other realizations of the initial perturbation ε(t0) for AR = 0.15 and Δx = 1 m. The linear estimate is plotted with symbols. Underprediction of the actual growth is robust.

Fig. 5.

Other realizations of the initial perturbation ε(t0) for AR = 0.15 and Δx = 1 m. The linear estimate is plotted with symbols. Underprediction of the actual growth is robust.

4. Spatial structure of the errors

In this section the applicability of linear theory to the spatial structure of the error fields is examined for AR = 0.15, 0.60, and 1.50.

Time-averaged relative errors 〈〉(x, z) are shown in Fig. 6. The structure of the error fields depends on AR. There is essentially complete loss of predictability inside the canopy for AR = 0.15. This is the behavior that is widely studied in analyses of homogeneous turbulence (Leith and Kraichnan 1972) and global atmospheric models (Tribbia and Baumhefner 2004; Ngan and Eperon 2012). Residual predictability remains for deeper canyons—for example, AR = 0.60 and 1.50; this case was extensively analyzed in Lo and Ngan (2015).

Fig. 6.

Time-averaged relative errors 〈〉(x, z) for AR = (a) 0.15, (b) 0.60, and (c) 1.50. Data were averaged over y and from τsat to the end of the simulation (yielding ~370–450 samples). Note that nondimensional coordinates are used; hence the appearance of the buildings does not change with AR. For clarity, the values are capped at 〈〉 = 2.

Fig. 6.

Time-averaged relative errors 〈〉(x, z) for AR = (a) 0.15, (b) 0.60, and (c) 1.50. Data were averaged over y and from τsat to the end of the simulation (yielding ~370–450 samples). Note that nondimensional coordinates are used; hence the appearance of the buildings does not change with AR. For clarity, the values are capped at 〈〉 = 2.

Time-averaged raw error fields 〈〉(x, z) are shown in Fig. 7. In all cases, the largest errors are found in the shear layer, where shear instability produces strong turbulence (e.g., Letzel et al. 2008). Differences among the three aspect ratios are easily discerned. For AR = 0.15, there is significant penetration of errors into the canyon; for AR = 0.60, errors extend well above the roof level; for AR = 1.50, errors are more tightly localized.

Fig. 7.

Time-averaged error fields 〈〉(x, z) (m2 s−2) for AR = (a) 0.15, (b) 0.60, and (c) 1.50. The averaging procedure is identical to that in Fig. 6. The structure of the raw error fields depends on the aspect ratio and flow regime.

Fig. 7.

Time-averaged error fields 〈〉(x, z) (m2 s−2) for AR = (a) 0.15, (b) 0.60, and (c) 1.50. The averaging procedure is identical to that in Fig. 6. The structure of the raw error fields depends on the aspect ratio and flow regime.

The interpretation of the relative errors can be misleading for inhomogeneous flows. Regions with relatively weak flow and small values of E (e.g., corners) may have large relative errors even though the raw error ΔE is small. It is debatable whether this effect translates to loss of predictability in a meaningful sense. Henceforth we focus on r.

During the linear growth stage (τ < τsat), it is plausible that the structure of the errors should reflect the largest eigenvalues. Errors should be largest in regions with strong strain, such as the shear layer. A striking resemblance exists between the time-averaged eigenvalues (Fig. 8) and errors (Fig. 7); 〈〉 is largest where 〈〉 is largest. Moreover, the same characteristic structures can be seen in both figures (e.g., shear layer and canyon vortex).

Fig. 8.

Time-averaged eigenvalues 〈〉(x, z) (s−1) for AR = (a) 0.15, (b) 0.60, and (c) 1.50. The averaging procedure is identical to that in Fig. 6. The spatial structure of 〈〉 resembles that of 〈〉.

Fig. 8.

Time-averaged eigenvalues 〈〉(x, z) (s−1) for AR = (a) 0.15, (b) 0.60, and (c) 1.50. The averaging procedure is identical to that in Fig. 6. The spatial structure of 〈〉 resembles that of 〈〉.

The strong resemblance between 〈〉 and 〈〉 suggests that the statistical structure of the error field after saturation is controlled by the linear dynamics. It is possible for linear dynamics to remain important even after saturation because urban flows are inhomogeneous and perturbations do not grow uniformly along Lagrangian trajectories. While perturbations grow rapidly near the roof level, where the strain is strongest, they do not stay confined to the shear layer; instead they may enter the canopy or be swept away from the canyon, where the perturbation decays. This is supported by the time-averaged error fields (Fig. 7), which indicate that the error decays rapidly away from the roof level. Animations of the error field (cf. Lo and Ngan 2015) indicate that errors do not continue to grow after fluid parcels leave the shear layer. Strong shear near the roof level is the driving force behind the error growth, but shear elsewhere may serve to inhibit growth. The role of shear in suppressing instabilities is well known in fluid dynamics (cf. Dritschel et al. 1991).

The foregoing may seem counterintuitive since the linear term does not dominate for τ > τsat. Note, however, that there are fundamental differences between street-canyon flows and homogeneous turbulence or large-scale atmospheric dynamics. In the latter case, strong nonlinearity throughout the domain generates an inverse cascade of errors (cf. Ngan et al. 2009). In the former case, nonlinearity dominates only in the shear layer (cf. Figs. 2 and 3). Nonlinearity does not destroy the error pattern established by the linear dynamics because the turbulence is strongly inhomogeneous.

Spatial correlations between 〈〉 and 〈〉 have been calculated using Eq. (8). As expected, the correlation between the time-averaged eigenvalues and error is strong: C(〈〉, 〈〉) = 0.83, 0.90, and 0.92 for AR = 0.15, 0.60, and 1.50 respectively.

Although the relative error is O(1) over a significant fraction of the canopy, linear theory is still relevant to the long-time error (and by extension, the long-time evolution of the perturbed flow) in a statistical sense. This is not possible in homogeneous turbulence because ensemble members become almost completely decorrelated on sufficiently long time scales.

Similar results are obtained when a continuous forcing is used in place of errors in the initial conditions. Figure 9 shows that 〈〉 is nearly identical for perturbations to the initial conditions and the time-dependent fields. This is consistent with the claim that the errors are determined by the linear flow. Time series of the relative error are virtually indistinguishable from Fig. 1 (not shown).

Fig. 9.

The 〈〉 (m2 s−3) for small errors in the (a) initial conditions and (b) time-dependent model fields. In the latter case, AR = 0.40, Δx = 1 m, and the perturbation is applied continuously. The time-averaged error field does not depend on the nature of the forcing.

Fig. 9.

The 〈〉 (m2 s−3) for small errors in the (a) initial conditions and (b) time-dependent model fields. In the latter case, AR = 0.40, Δx = 1 m, and the perturbation is applied continuously. The time-averaged error field does not depend on the nature of the forcing.

5. Aspect-ratio dependence

The previous sections have highlighted the behavior for aspect ratios that are representative of the isolated-roughness, wake-interference, and skimming-flow regimes. We now consider how time-averaged correlations with ΔE depend on AR.

Figure 10a shows correlations between the time-dependent error and the time-averaged eigenvalues: . For both medium (Δx = 1m) and high (Δx = 0.5 m) resolution, C is greater than approximately 0.75. This implies that the time-averaged eigenvalues carry significant information about the error field after saturation, in agreement with the idea that nonlinearity does not destroy the statistical structure of the error field (established prior to saturation).

Fig. 10.

Time-averaged correlation between the raw error and the largest eigenvalues: (a) time-averaged eigenvalues and (b) instantaneous eigenvalues . The labels in the former case correspond to medium (mr; Δx = 1 m) and high (hr; Δx = 0.5 m) resolution; the lines in the latter case correspond to different snapshots (i.e., choices of τ) for low resolution. The ΔE was calculated for τ > τsat, whereas 〈〉 represents an average over τ > 0.

Fig. 10.

Time-averaged correlation between the raw error and the largest eigenvalues: (a) time-averaged eigenvalues and (b) instantaneous eigenvalues . The labels in the former case correspond to medium (mr; Δx = 1 m) and high (hr; Δx = 0.5 m) resolution; the lines in the latter case correspond to different snapshots (i.e., choices of τ) for low resolution. The ΔE was calculated for τ > τsat, whereas 〈〉 represents an average over τ > 0.

The ΔE is also correlated with the instantaneous eigenvalues. Figure 10b shows , correlations between ΔE(τ′) and the eigenvalues evaluated at specific times λ(τ′). Although the values are generally smaller than those for the time-averaged eigenvalues , the correlation is still fairly strong, with C being greater than approximately 0.6. The associated p values, which indicate the probability that a similar correlation could occur by chance, are extremely small; the spatial correlations represent an average over ~105 points. The aspect-ratio dependence is stronger than in Fig. 10a: there is a local maximum around AR = 0.40 and a local minimum around AR = 0.65. The trend is similar for Δx = 0.5 m (not shown).

The instantaneous correlations are weaker. Because the control velocity Ui is turbulent, the approximation to the velocity gradient

 
formula

is inherently noisy and may not be well correlated with the time average

 
formula

everything else being the same, the correlation between ∂Ui/∂xj and should increase when Ui is more spatially coherent. Hence differences between and should depend on aspect ratio and flow regime.

This claim is supported by Fig. 10b. The local maximum occurs at an aspect ratio that is close to the transition value from isolated roughness to wake interference for an isolated street canyon (AR of ~0.3–0.4); the local minimum occurs close to the transition value from wake interference to skimming flow (AR of ~0.65–0.75). At the critical transition values, long-range interactions change character (Ngan and Lo 2016), and it is natural that this be reflected in ∂Ui/∂xj and . Indeed, the temporal variability

 
formula

where Es and E denote spatial and temporal averages, is largest for AR = 0.60 (not shown). The differences are relatively small, however (about 20% when compared with AR = 0.15 and 1.50). For a deep street canyon with AR > 1, the correlations are strong (C greater than approximately 0.6) and linear theory is most relevant. This is consistent with the idea that turbulence within a deep street canyon weakens as the aspect ratio increases (cf. Li et al. 2006).

It is possible that correlation between the error and eigenvalue fields could be a by-product of a strong correlation between ΔE and U. For AR of less than approximately 1, however, the correlation is weak; has a maximum value of 0.24 (not shown). This suggests that the correlation between ΔE and λ is not fortuitous and that the (linear) strain field has more influence on the error evolution than the (linear) streamwise velocity. For AR = 1.50, however, the difference between and is much smaller. Nevertheless, both correlations are large.

6. Discussion

Errors and error growth inevitably accompany any forecast. Hence a practical operational system requires a means to correct for them. The feasibility and ultimate success of such a system depend on an appropriate model of error growth. In operational data assimilation, a tangent linear model (TLM) is commonly adopted (Courtier et al. 1994).

This work demonstrates that a linear model error can also be applied to turbulent flow in a street canyon. First, the error growth rate can be estimated, with errors around 10%–40%, from the largest eigenvalues of the velocity gradient of the unperturbed linear flow. This supports the validity of the linear perturbation Eq. (10), which is not strictly identical to a TLM. Because the linearized pressure term is not retained, it represents a linear “perturbation forecast” model in which knowledge of the physics is exploited (Rawlins et al. 2007). Second, the linear error model is relevant even after nonlinearity becomes important and the errors saturate. This is evidenced by the strong spatial correlations between the error and eigenvalue fields, time averaged and instantaneous: some predictability is retained in a statistical sense. Such behavior is very different from that of large-scale atmospheric dynamics and homogeneous turbulence (cf. Tanguay et al. 1995), wherein strong nonlinearity throughout the domain causes data-assimilation methods that are based on a TLM to fail on time scales that are comparable to or greater than the nonlinear time scale. In turbulent street-canyon flow, the growth mechanism is essentially advective as the strongest growth occurs when the perturbation is aligned with the strain; moreover, the presence of an open boundary limits growth.

The practical implication is that standard data-assimilation algorithms that are based on linearization, such as the extended Kalman filter and incremental 4D-Var, should work well for turbulent flow inside street canyons. This would accelerate development and implementation. Nonlinear ensemble data-assimilation methods have theoretical advantages (Lorenc 2003), but the computational cost of a building-resolving ensemble would be extremely high. If operational building-resolving forecasts were to become a reality, less computationally demanding approaches—for example, linear data-assimilation algorithms—would almost certainly be required in the first instance.

A related implication is that a linear error or perturbation model could be used to improve the modeling of urban canopy winds. In the absence of time-dependent boundary conditions or forcing, the ability of a building-resolving model to forecast the time evolution of urban flow is limited. A linear perturbation model could be used to represent the effects of rapid variations (in, e.g., inflow winds) on the time-averaged state. This would be useful for applications, such as dispersion modeling, in which coupling and operational cycling are not feasible. The best-known theoretical models of urban canopy winds also assume steady boundary conditions (Cionco 1965; Macdonald 2000; Di Sabatino et al. 2008).

The linear error model does not apply equally well to all urban geometries. The linear estimate of the growth rate (Fig. 1) and the spatial correlation between eigenvalues and errors (Fig. 10) depend on the aspect ratio. With respect to , agreement is better for the wake-interference and skimming-flow regimes. With respect to , the correlation is greater in the isolated-roughness and skimming-flow regimes. The upshot is that linear theory is most relevant to skimming flow and deep street canyons with AR of greater than approximately 1. In this regime, turbulence within the canopy weakens as the separation increases between the shear layer, where the turbulence is generated, and the bottom of the canopy.

The direct applicability of this paper to real urban areas rests on several assumptions. First, the canopies must be sufficiently deep for linear theory to be accurate. Second, real canyons can be modeled as idealized 2D canyons without intersections. Third, thermal effects can be neglected.

The first two assumptions concern the geometry of real urban canyons. Deep street canyons may be found in areas with many skyscrapers, for example, Manhattan in New York City, New York, and Central in Hong Kong, China. Canyons with AR of greater than approximately 1 are also found in neighborhoods that are densely packed although not necessarily that have the tallest buildings. A wider range of urban geometries and cities could be considered if lower accuracy were sufficient. Calculations with a street-canyon array in which the length-to-separation ratio of the canyons is 8:1 indicate that the spanwise-averaged flow resembles that of a single street canyon away from the intersection (not shown). Wind-tunnel results indicate that the sensitivity to the canyon length is fairly weak for deep canyons in the skimming-flow regime and with ratios of 2:1 and greater (Oke 1988).

As a concrete example, we consider a small area within Mong Kok, a mixed-use Hong Kong neighborhood with commercial and residential buildings. Analyses of digital building data for this area show that smaller side streets have width W of approximately 10–15 m, lengths L of approximately 50 m, and building heights H of approximately 40 m. The main thoroughfare (Nathan Road) has W of approximately 30 m, H of approximately 42 m, L of approximately 60 m, and an average aspect ratio of ~1.45 (Lo and Ngan 2017). Assuming the height of a story may be taken to be 3 m, buildings in this area are around 13 stories, which is unremarkable by Hong Kong standards but rather tall for many European cities. Many urban areas in Asia and North America contain buildings of this height and taller. If one assumes that the canyon width and separation are the same, then L/s is approximately 3–5 and 1.5 for the side streets and thoroughfare, respectively. This suggests that real urban canyons may not be sufficiently long for the results of this paper to be directly applicable.

The second assumption means that error growth resulting from the heating of buildings or boundaries is excluded. It seems unlikely that the error growth will be driven purely by advection and strain when there is strong heating.

The preceding assumptions limit the applicability of the results described in this paper. Nevertheless, it may be possible to extend them if a regime analogous to skimming flow in a 2D canyon is generalized to realistic geometry and thermal forcing. It is necessary that the flow be spatially coherent and approximately steady (in the Lagrangian frame). Extending the eigenvalue problem Eq. (10) to a general linear stability analysis is straightforward, but further work will be needed required to verify the underlying assumptions. Despite its limitations, the work presented here applies to an important limit (deep canyons and neutral flow) that is still being actively researched. It should represent a useful reference point for future studies.

Acknowledgments

This research was supported by the Research Grants Council of Hong Kong (Project 21304515) and City University of Hong Kong (Projects 7004165 and 7200403). We are grateful to the anonymous reviewers for many valuable comments and suggestions. Nicole Kong assisted with the analyses of the building data.

REFERENCES

REFERENCES
Barlow
,
J.
,
2014
:
Progress in observing and modelling the urban boundary layer
.
Urban Climate
,
10
,
216
240
, doi:.
Cionco
,
R. M.
,
1965
:
A mathematical model for air flow in a vegetative canopy
.
J. Appl. Meteor.
,
4
,
517
522
, doi:.
Clayton
,
A. M.
,
A. C.
Lorenc
, and
D. M.
Barker
,
2013
:
Operational implementation of a hybrid ensemble/4D-Var global data assimilation system at the Met Office
.
Quart. J. Roy. Meteor. Soc.
,
139
,
1445
1461
, doi:.
Courtier
,
P.
,
J.-N.
Thépaut
, and
A.
Hollingsworth
,
1994
:
A strategy for operational implementation of 4DVar, using an incremental approach
.
Quart. J. Roy. Meteor. Soc.
,
120
,
1367
1387
, doi:.
Cui
,
Z.
,
X.
Cai
, and
C. J.
Baker
,
2004
:
Large-eddy simulation of turbulent flow in a street canyon
.
Quart. J. Roy. Meteor. Soc.
,
130
,
1373
1394
, doi:.
Daley
,
R.
,
1993
: Atmospheric Data Analysis. Cambridge University Press, 457 pp.
Deardorff
,
J. W.
,
1980
:
Stratocumulus-capped mixed layers derived from a three-dimensional model
.
Bound.-Layer Meteor.
,
18
,
495
527
, doi:.
Di Sabatino
,
S.
,
E.
Solazzo
,
P.
Paradisi
, and
R.
Britter
,
2008
:
A simple model for spatially-averaged wind profiles within and above an urban canopy
.
Bound.-Layer Meteor.
,
127
,
131
151
, doi:.
Dritschel
,
D. G.
,
P. H.
Haynes
,
M. N.
Juckes
, and
T. G.
Shepherd
,
1991
:
Stability of a two-dimensional vorticity filament under uniform strain
.
J. Fluid Mech.
,
230
,
647
665
, doi:.
Durran
,
D. R.
, and
M.
Gingrich
,
2014
:
Atmospheric predictability: Why butterflies are not of practical importance
.
J. Atmos. Sci.
,
71
,
2476
2488
, doi:.
Foken
,
T.
,
2008
: Micrometeorology. Springer-Verlag, 306 pp., doi:.
Franke
,
J.
,
A.
Hellsten
,
H.
Schlünzen
, and
B.
Carissimo
, Eds.,
2007
: Best practice guideline for the CFD simulation of flows in the urban environment—COST Action 732: Quality assurance and improvement of microscale meteorological models. COST Office Tech. Rep., 52 pp. [Available online at http://theairshed.com/pdf/COST%20732%20Best%20Practice%20Guideline%20May%202007.pdf.]
Hanna
,
S. R.
, and Coauthors
,
2006
:
Detailed simulations of atmospheric flow and dispersion in downtown Manhattan: An applications of five computational fluid dynamics models
.
Bull. Amer. Meteor. Soc.
,
87
,
1713
1726
, doi:.
Hunter
,
L.
,
I.
Watson
, and
G.
Johnson
,
1990
:
Modelling air flow regimes in urban canyons
.
Energy Build.
,
15
,
315
324
, doi:.
Inagaki
,
A.
,
M.
Castillo
,
Y.
Yamashita
,
M.
Kanda
, and
H.
Takimoto
,
2012
:
Large-eddy simulation of coherent flow structures within a cubical canopy
.
Bound.-Layer Meteor.
,
142
,
207
222
, doi:.
Kalnay
,
E.
,
2002
: Atmospheric Modeling, Data Assimilation and Predictability. Cambridge University Press, 341 pp.
Leith
,
C.
, and
R. H.
Kraichnan
,
1972
:
Predictability of turbulent flows
.
J. Atmos. Sci.
,
29
,
1041
1057
, doi:.
Letzel
,
M. O.
,
M.
Krane
, and
S.
Raasch
,
2008
:
High resolution urban large-eddy simulation studies from street canyon to neighbourhood scale
.
Atmos. Environ.
,
42
,
8770
8784
, doi:.
Letzel
,
M. O.
,
C.
Helmke
,
E.
Ng
,
X.
An
,
A.
Lai
, and
S.
Raasch
,
2012
:
LES case study on pedestrian level ventilation in two neighbourhoods in Hong Kong
.
Meteor. Z.
,
21
,
575
589
, doi:.
Li
,
X.-X.
,
C.-H.
Liu
,
D. Y. C.
Leung
, and
K. M.
Lam
,
2006
:
Recent progress in CFD modelling of wind field and pollutant transport in street canyons
.
Atmos. Environ.
,
40
,
5640
5658
, doi:.
Liu
,
C.-H.
, and
M.
Barth
,
2002
:
Large-eddy simulation of flow and scalar transport in a modeled street canyon
.
J. Appl. Meteor.
,
41
,
660
673
, doi:.
Lo
,
K. W.
, and
K.
Ngan
,
2015
:
Predictability of turbulent flow in street canyons
.
Bound.-Layer Meteor.
,
156
,
191
210
, doi:.
Lo
,
K. W.
, and
K.
Ngan
,
2017
:
Characterizing ventilation and exposure in street canyons using Lagrangian particles
.
J. Appl. Meteor. Climatol.
,
56
,
1177
1194
, doi:.
Lorenc
,
A. C.
,
2003
:
The potential of the ensemble Kalman filter for NWP—A comparison with 4D-Var
.
Quart. J. Roy. Meteor. Soc.
,
129
,
3183
3203
, doi:.
Lorenz
,
E. N.
,
1969
:
The predictability of a flow which possesses many scales of motion
.
Tellus
,
21
,
289
307
, doi:.
Macdonald
,
R. W.
,
2000
:
Modelling the mean velocity profile in the urban canopy layer
.
Bound.-Layer Meteor.
,
97
,
25
45
, doi:.
Maronga
,
B.
, and Coauthors
,
2015
:
The Parallelized Large-Eddy Simulation Model (PALM) version 4.0 for atmospheric and oceanic flows: Model formulation, recent developments, and future perspectives
.
Geosci. Model Dev.
,
8
,
2515
2551
, doi:.
Ngan
,
K.
, and
G. E.
Eperon
,
2012
:
Middle atmosphere predictability in a numerical weather prediction model: Revisiting the inverse error cascade
.
Quart. J. Roy. Meteor. Soc.
,
138
,
1366
1378
, doi:.
Ngan
,
K.
, and
K. W.
Lo
,
2016
:
Revisiting the flow regimes for urban street canyons using the numerical Green’s function
.
Environ. Fluid Mech.
,
16
,
313
334
, doi:.
Ngan
,
K.
,
D. N.
Straub
, and
P.
Bartello
,
2004
:
Three-dimensionalization of freely-decaying two-dimensional turbulence
.
Phys. Fluids
,
16
,
2918
2932
, doi:.
Ngan
,
K.
,
P.
Bartello
, and
D. N.
Straub
,
2009
:
Predictability of rotating stratified turbulence
.
J. Atmos. Sci.
,
66
,
1384
1400
, doi:.
Oke
,
T.
,
1988
:
Street design and urban canopy layer climate
.
Energy Build.
,
11
,
103
113
, doi:.
Palmer
,
T.
, and
R.
Hagedorn
, Eds.,
2006
: Predictability of Weather and Climate. Cambridge University Press, 702 pp.
Park
,
S.-B.
, and
J.-J.
Baik
,
2013
:
A large-eddy simulation study of thermal effects on turbulence coherent structures in and above a building array
.
J. Appl. Meteor. Climatol.
,
52
,
1348
1365
, doi:.
Raasch
,
S.
, and
M.
Schröter
,
2001
:
PALM—A large-eddy simulation model performing on massively parallel computers
.
Meteor. Z.
,
10
,
363
372
, doi:.
Rawlins
,
F.
,
S. P.
Ballard
,
K. J.
Bovis
,
A. M.
Clayton
,
D.
Li
,
G. W.
Inverarity
,
A. C.
Lorenc
, and
T. J.
Payne
,
2007
:
The Met Office global four-dimensional variational data assimilation scheme
.
Quart. J. Roy. Meteor. Soc.
,
133
,
347
362
, doi:.
Schlunzen
,
K.
,
D.
Grawe
,
S.
Bohnenstengel
,
I.
Schluter
, and
R.
Koppmann
,
2011
:
Joint modelling of obstacle induced and mesoscale changes—Current limits and challenges
.
J. Wind Eng. Ind. Aerodyn.
,
99
,
217
225
, doi:.
Stauffer
,
D. R.
, and
N. L.
Seaman
,
1994
:
Multiscale four-dimensional data assimilation
.
J. Appl. Meteor.
,
33
,
416
434
, doi:.
Tanguay
,
M.
,
P.
Bartello
, and
P.
Gauthier
,
1995
:
Four-dimensional data assimilation with a wide range of scales
.
Tellus
,
47A
,
974
997
, doi:.
Tribbia
,
J. J.
, and
D. P.
Baumhefner
,
2004
:
Scale interactions and atmospheric predictability: An updated perspective
.
Mon. Wea. Rev.
,
132
,
703
713
, doi:.
Yamada
,
T.
, and
K.
Koike
,
2011
:
Downscaling mesoscale meteorological models for computational wind engineering applications
.
J. Wind Eng. Ind. Aerodyn.
,
99
,
199
216
, doi:.
Zajaczkowski
,
F. J.
,
S. E.
Haupt
, and
K. J.
Schmehl
,
2011
:
A preliminary study of assimilating numerical weather prediction data into computational fluid dynamics models for wind prediction
.
J. Wind Eng. Ind. Aerodyn.
,
99
,
320
329
, doi:.

Footnotes

© 2017 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).

1

To be more precise, long-range interactions between the roof-level and pedestrian-level shear layers take extremal values.

2

An empirical criterion that was based on the second derivative was tested, but it requires manual adjustment of thresholds.