# Laboratory Simulations of Local Winds in the Atmospheric Boundary Layer via Image Analysis

Authors: Monica Moroni; Antonio Cenedese
Journal: Advances in Meteorology
 (2015)

---

## Abstract

In the atmospheric boundary layer, under high pressure conditions and negligible geostrophic winds, problems associated with pollution are the most critical. In this situation local winds play a major role in the evaluation of the atmospheric dynamics at small scales and in dispersion processes. These winds originate as a result of nonuniform heating of the soil, either when it is homogeneous or in discontinuous terrain in the presence of sea and/or slopes. Depending on the source of the thermal gradient, local winds are classified into convective boundary layer, sea and land breezes, urban heat islands, and slope currents. Local winds have been analyzed by (i) simple analytical models; (ii) numerical models; (iii) field measurements; (iv) laboratory measurements through which it is impossible to completely create the necessary similarities, but the parameters that determine the phenomenon can be controlled and each single wind can be separately analyzed. The present paper presents a summary of laboratory simulations of local winds neglecting synoptic winds and the effects of Coriolis force. Image analysis techniques appear suitable to fully describe both the individual phenomenon and the superposition of more than one local wind. Results do agree with other laboratory studies and numerical experiments.


---

## Body

## 1. Introduction

The atmospheric boundary layer (ABL) is the portion of the atmosphere directly influenced by the earth’s surface with a characteristic response time to forcing on the order of less than one hour [1, 2]. It is the region closer to the ground where thermodynamic parameters depend on the daily evolution of the solar radiation. The ABL thickness may vary from a few hundred to 2000 meters. Many human activities occur in the ABL and the pollutants produced will remain for long periods of times. An investigation of the ABL is therefore mandatory both for weather prediction (the weather is strictly related to temperature, humidity, and wind) and for air quality issues.Even when no geostrophic winds occur in the ABL, air movements (hereinafter called local winds [3]) may be generated due to the baroclinic behavior of the fluid [4, 5]. Then, the pressure gradient ∇p is linked to the temperature and density gradients, ∇T and ∇ρ. Depending on the origin of the temperature gradient triggering the air movement, local winds are classified as follows.Convective Boundary Layer. This develops as a direct consequence of solar radiation during fair weather conditions and negligible geostrophic winds; mixed layers are usually daytime phenomena, with growth caused by entrainment of free-atmosphere air into the mixed-layer top; updraft and downdraft structures with horizontal and vertical velocity components of the same order of magnitude characterize the air movement.Sea and Land Breeze. The large heat capacity of water is responsible for the constancy of the sea temperature during the diurnal cycle whereas the ground temperature may vary by more than ten degrees in the same time frame, generating a horizontal temperature gradient and a corresponding parallel wind, oriented from the sea toward the land during the day and from the land toward the sea during the night; compensating currents take place higher up in the atmosphere; similar circulations can establish close to large lakes or rivers as well.Urban Heat Island. Due to the low albedo, urban areas tend to have a temperature some degrees higher than neighboring rural zones; this triggers air movement toward the city close to the ground and a return current at larger heights.Slope Currents. Due to the horizontal temperature difference between air adjacent to a slope and the ambient air at the same altitude, the flow is upslope (or anabatic) during the daytime when the slope is warmer and downslope (or katabatic) during the nighttime when the slope temperature drops; if a valley is included within the slopes and differential heating conditions occur, cross-valley and/or along-valley wind systems, perpendicular and longitudinal to the valley axis, respectively, may occur; air motions close to the ground are compensated by a return flow higher up in the atmosphere.It is, in general, unlikely to isolate each local wind because they usually interact with each other and with geostrophic winds.Local winds are typically unsteady phenomena. The quantities influencing and/or describing air movement can be divided into a mean part (ensemble average) and a perturbation part related to either a wave effect and turbulence effect. The phenomenon is usually assumed to be instantaneously ergodic and ensemble averages are evaluated as temporal means over time intervals of some order of magnitude shorter than the phenomenon time scale. This implies applying a low-pass filter which allows isolation of large scale vortices from turbulence. The presence of waves can be found through Fourier and/or correlation analysis capable of recognizing space and time periodicity [2].Local winds have been studied by analytical modeling, numerical simulations, laboratory experiments, and field surveys. None of these can provide completely exhaustive answers. Therefore their combined use appears to be the most appropriate way to characterize such a complex phenomenon.Analytical solutions and bulk models [6–9] use radical simplifications to provide a description, sometimes only qualitative, of the phenomenon, elucidating its most relevant physical aspects.The numerical modeling of turbulent flows has been widely employed due to the ever growing capacity of computers. However, a direct simulation (direct numerical simulation, DNS) capable of describing the phenomenon at all scales appears impractical since the smallest allowable cell dimensions would be still much larger than the smallest turbulent structures. The models used are generally based on the Reynolds decomposition in mean quantities and fluctuating components (Reynolds-averaged Navier-Stokes, RANS) parameterizing the turbulence closure through empirical or semiempirical relations. Limited area models (LAMs) are used to study local winds coupled to general circulation models to obtain the necessary boundary conditions. LAMs provide results with a resolution on the order of kilometers by adopting the nested cell technique [10]. More recently, higher resolutions (on the order of hundreds of meters) can be achieved through the LES (large eddy simulation) technique which filters turbulence at different scales [11]. An alternative method to study even smaller scales is to use CFD techniques, which describe only limited spatial portions. As highlighted by Wyngaard [12], there is an area still not fully explored, “Terra Incognita.” This area is of great interest for the study of local winds, where the parameterization of the turbulence is not completely satisfactory.The field measurement campaigns, using probes placed on masts, towers, aircraft, and weather balloons or through remote sensing [13, 14], are the only way to get a real description of the phenomenon. However, these measurements are particularly onerous and provide a field description in a limited number of points only.Laboratory simulations allow control of the causes of local winds and thus present the possibility of analyzing them separately. High quality and detailed datasets for a range of convective conditions are mandatory as test cases for checking the increasingly sophisticated numerical models [15]. The most significant limitation of these simulations is the inability to achieve all similarity conditions, since Reynolds and/or Rayleigh numbers in the laboratory simulations are a few orders of magnitude lower than the real values. However, for sufficiently large values of these parameters, it can be assumed that the phenomenon does not depend on the parameters themselves.Given the difficulty in their use, only a limited number of wind tunnels have been used for the study of the atmospheric boundary layer in not neutral conditions [16–21]. The implementation of thermal stratification is particularly difficult and expensive in terms of energetic demand. Furthermore, it is difficult, if not impossible, to simulate strongly convective conditions (with low tunnel speeds) and maintain suitable scaling [22]. The advantages of wind tunnels reside in their ability to model the effects of local topography as well as the combined effects of wind shear and convective turbulence at high wind speeds, thus complementing the zero to light wind regime that is best simulated in water tank models.The experiments are more numerous in water tanks where the time evolution of the convective boundary layer (light winds) can be modeled. The initial fluid stratification can be achieved through different salt concentrations [23] or different temperatures [24, 25]. Among the advantages of saline convection over thermal convection is the absence of heat loss through the walls of the tank. On the other hand, Hibberd and Sawford [15, 23] describe difficulties in controlling the buoyancy flux and obtaining repeatable results as well as ensuring a constant flux over the ground surface. Experiments do allow one to fully characterize the fluid dynamic field through visualizations, temperature, and multipoint velocity measurements.In this work, we will refer to investigations developed in the Hydraulic Laboratory of DICEA-Sapienza University of Rome. Experiments have been conducted in test sections filled with distilled water. Distilled water is used as working fluid, to allow both a large heating rate and sufficient time to take measurements of the evolving thermal structures. The tank sides are made of glass to ensure optical access. A polystyrene sheet is placed over the top of the tank to reduce heat loss. The temperature of the sidewalls is not controlled but can be reasonably considered constant throughout the experiment. Pollen particles of about 80μm mean diameter, assuming the same density as the seeded liquid [26–28], are used as neutrally buoyant passive tracer to reconstruct particle trajectories. The test section is illuminated through a planar light sheet obtained through a laser or a high power lamp. Images are acquired with 1 Mbyte, 8-bit pixels camera at 25 Hz; due to the low velocity, images have been occasionally undersampled. The velocity field has been determined through Lagrangian particle tracking techniques, which allow reconstructing particle trajectories from which velocity and acceleration can be obtained directly. The sparse velocity vectors are interpolated on a regular grid to gather an Eulerian description of the flow field. Temperature is detected through T-type (cupper-constantan) thermocouples of uncertainty less than 0.1 K and sample frequency of 1 Hz. Thermocouples are placed within the test section along a vertical array to measure vertical profiles and on the lower boundary to test horizontal homogeneity in supplying heat. In both cases thermocouples do not disturb the flow field since their position is few centimeters far away from the illuminated plane. Heat exchanges at the boundary surfaces are thermocontrolled through thermocryostats, Peltier cells, and electric heaters. If available, the laboratory results have been compared to numerical simulation and field measurements.The following assumptions apply:(i)
Coriolis forces being negligible;(ii)
geostrophic winds being absent;(iii)
incompressible fluid and validity of the Boussinesq hypothesis;(iv)
2D flow field (in general planar motion), axial-symmetry for the heat islands;(v)
no humidity effects: that is, the virtual potential temperatureθv is the same as the potential temperature θ, which is in turn the temperature T when experiments are conducted in water.The notation for all the cases under investigation has been kept consistent when possible. Quantities are defined the first time they appear in the text whereas the most common quantities are in the list of symbols.The paper is organized as follows: Section2 describes image analysis techniques; Section 3 focuses on the description of the experimental facilities and main results of the convective boundary layer; Section 4 describes the urban heat island experiments and main results; Section 5 presents the laboratory model to reproduce anabatic and katabatic currents and the corresponding analytical models; Section 6 presents one case of interaction between local winds, slope flows, and urban heat islands; the concluding section draws the main outcomes of the experimental efforts and comparisons.


## 2. Flow Field in a Fluid Measured through Image Analysis Techniques

The need for measuring velocity fields has historically led researchers to develop experimental techniques and related instrumentation. An “ideal” measurement system should be nonintrusive to avoid flow field perturbations, should not require calibration, and should be suitable for obtaining the velocity field with a time and space resolution smaller than the characteristic time and length scales (i.e., Kolmogorov scale for turbulent flows). Particle image velocimetry (PIV) and Lagrangian particle tracking (LPT), based on optical methods capable of providing the velocity of image regions, represent the best approximation of this “ideal” system. These techniques allow simultaneous multipoint velocity measurements in a plane (2D measurement techniques) or in a volume (3D measurement techniques). PIV evaluates the average displacement of particles belonging to a subregion (interrogation area) of the whole imaged area using correlation between couple of images. PIV reconstructs the Eulerian velocity field on a regular, equispaced grid. LPT includes all the methods in which the displacements of every single particle are considered to reconstruct trajectories. Velocity is then evaluated from a Lagrangian point of view considering the ratio between particle displacement along the reconstructed trajectories and the correspondent time interval. LPT provides sparse velocity vectors at points coincident with particle centroid positions. In this work, different 2D LPT techniques have been used, namely, feature tracking (FT [26]) and hybrid Lagrangian particle tracking (HLPT [29]).The Lagrangian description is particularly suitable to study dispersion phenomena. The transition from one description to another is always possible. While shifting from the Lagrangian to the Eulerian description implies an integration, which allows the use of a low pass filter to reduce noise, the transition from the Eulerian to the Lagrangian description requires the solution of a system of differential equations that amplifies noise.Feature tracking and hybrid Lagrangian particle tracking share the principal steps required to analyze images acquired at fixed time intervals. As an example of the measurement chain, we will examine an anabatic wind experiment:(i)
acquiring images by using standard cameras (25 frames/s, 576 × 764 pixels) or higher spatial resolution (up to 4 megapixels) and/or time resolution (up to 1000 frames/s) cameras;(ii)
noise attenuation and elimination of boundaries and areas with no particles;(iii)
recognition of a particle and the coordinates of its centroid.The particle detection algorithm used to analyze the images acquired during the experiments reported herein is based on the optical flow equation, which defines the conservation of the image intensity [29];(iv)
particle tracking and Lagrangian description of the flow field.Various tracking techniques have been developed to identify successive positions of the same particle and thus extract the displacement and velocity of that particle along its trajectory [30] (Figure 1);(v)
reconstruction of the Eulerian flow field.Postprocessing and presentation of experimental data are employed to map randomly spaced vector fields onto regular grids, to replace erroneous vectors with values computed from the neighboring vectors or to refine the original grid (Figure2);(vi)
evaluation of spatial derivatives, divergence, and vorticity (Figure2);(vii)
streamlines reconstruction.The analysis of streamlines provides a topological description of fluid-dynamic field highlighting the presence of critical points: nodes, saddles, and foci (Figure3).Figure 1
Particle centroids and trajectories reconstructed by HLPT: (a) whole image and (b) zoom in the area highlighted.
(a)
(b)Figure 2
Eulerian velocity fields and vorticity map (the vorticity unit is s−1).Figure 3
Streamlines.


## 3. Convective Boundary Layer

High pressure conditions and no clouds or geostrophic winds allow a stable layer to establish during the night due to the negative flux of sensible heat at the ground. Starting from the first hours of the morning, the flux of sensible heat at the ground becomes positive triggering an atmospheric unstable layer characterized by convective motions with the vertical velocity component of the same order of magnitude as the horizontal one. The fluid within this region, the height of which increases with time, is completely mixed.Taking into account the trend of the mean temperatureθ- and turbulent kinematic heat flux, define the following (Figure 4):(i)
zs height close to the ground where the mean potential temperature gradient changes from negative to zero;(ii)
ze height where the turbulent kinematic heat flux is zero;(iii)
zi inversion height where the turbulent kinematic heat flux is minimum; its value is 0.2 times the quantities at the lower boundary, w′θ′¯z=zi=-Aw′θ′¯z=0 with A=0.2;(iv)
zf minimum height where the turbulent kinematic heat flux and its derivatives are zero.Figure 4
Schematic representation of the temperature evolution within the atmospheric boundary layer at different time instants;θa indicates the background temperature; in (a) the superficial and entrainment layers are shown while they are not considered in (b).Different definitions have been used to define the same heights; for instance, the inversion heightzi can be defined as that where the second derivative of the mean temperature goes to zero (∂2θ-/∂z2=0 [31]) or turbulent intensities (u′2¯ or w′2¯) are a given fraction of the maximum value [31, 32]. The various definitions, however, lead to almost identical results.In the convective boundary layer, different regions can then be identified:(i)
the surface layer (SL) near the heat source has a negative temperature gradient meaning that an unstable condition acts as a driving force for the onset of convection; here temperature decreases fairly rapidly with height and heat transfer occurs mainly through conduction;(ii)
the mixing layer (ML), the bulk of the CBL where the mean vertical temperature gradients are nearly zero and the fluid is well mixed due to turbulence;(iii)
the entrainment zone (EZ), also referred to as the inversion layer or interfacial layer [33], a density interface between the turbulent (CBL) and nonturbulent region, where temperature increases significantly with height;(iv)
the free atmosphere (FA) or stable layer, which is not noticeably affected by the growing CBL and the fluid maintains its initial stratification.Indeed,w′u′¯ and w′θ′¯ go to zero at the surface in a layer (microlayer) with a thickness on the order of magnitude of the roughness due to vegetation and urban canopy or, for smooth boundaries, a thickness depending on where molecular diffusion dominates turbulent diffusion; in the microlayer heat transport is described by Fourier’s law.The vertical development of the CBL over homogeneous and flat continental areas has been described using tall masts and tethered balloons [13], aircraft [34, 35], or wind profilers [36, 37].Numerous water tank experiments have modeled the time evolution of the convective boundary layer. Field experiments aimed at measuring the turbulence budget of the CBL have shown that the mechanical generation of kinetic energy by wind shear is often confined close to the heat source supporting the validity of laboratory models in which no wind is present [38–40]. Willis and Deardorff [25, 41–43] and Deardorff and Willis [32, 44] systemized understanding of the convective boundary layer evolution and dispersion processes within it through laboratory experiments conducted in a heated water tank model (1.14 × 1.22 × 0.25 m3), initial fluid stratification achieved with different temperatures, measurements of horizontal and vertical temperature profiles through thermocouples, visual observations of penetrating thermals using a spread laser beam, and velocity measurements by use of streak photography of nonbuoyant particles. Kumar and Adrian [45] studied the entrainment zone and turbulence in a convection tank of dimensions 1.45 × 1.50 × 0.20 m3. Horizontally averaged temperature was measured with a resistance wire thermometer and thermocouples whereas laser Doppler velocimetry output velocity data. Hibberd and Sawford [15] present a saline rather than thermal convection. The 3.2 × 1.6 × 0.25 m3 tank has provided data for studies of turbulence, entrainment, and dispersion in the CBL via measurements of density by vertically-sounding single-point conductive probes and velocity by particle tracking velocimetry. Michaelian et al. [46] present experiments aimed at exploring the coupling of penetrative convection with internal waves in the adjoining stable layer. Tests were performed in a temperature-controlled, stably stratified fluid filling a long convection tank (0.845 × 0.077 × 0.107 m3). Temperature data were acquired by using temperature probe arrays whereas correlation image velocimetry was used to determine the velocity field associated with the flow.We will next describe the results obtained with a laboratory model designed to reproduce the penetrative shear-free convection phenomenon, in order to predict the CBL growth as a function of initial and boundary conditions, compute the dimensions of turbulent structures, understand the interaction between the turbulent and nonturbulent regions, and describe the fate of a passive tracer dissolved within the fluid through a transilient turbulence approach applied to a Lagrangian framework. Many experiments were performed [47, 48]. We will present next the results for a subset of two experiments (Table 1).Table 1
Convective boundary layer experimental parameters. Details of initial conditions for experiments A1 and A2.Tb0 is the temperature at the bottom before heating starts, TbC is the final heating temperature, and Γ is the vertical component of the temperature gradient. The duration of the experiments is roughly 2100 s.
Experiment number
T
b
0 (K)
T
b
C (K)
Γ
=
(
∂
θ
a
/
∂
z
) (K/m)
A1
283
308
101
A2
295
308
55


### 3.1. Experimental Procedure

The laboratory model consists of a convection chamber (0.40 × 0.40 × 0.41 m3) containing an initially stable, density stratified fluid, which is heated from below to trigger destabilization and penetrative convection. The fluid stratification is achieved through the two-tank methods [49]. The fluid initial conditions are velocity equal to zero and linear increase of temperature with height through a slope Γ. After the fluid stratification within the tank is completed, a thermostatically controlled hot water bath (set to a temperature Tc) is connected to the bottom metal plate and the experiment may begin. The framed area is Δx=0.200 m long (x-axis) and Δz=0.15 m high (z-axis).Thermocouples are placed within the test section along a vertical array of 26 probes to measure vertical profiles and on the lower boundary to test horizontal homogeneity in supplying heat. Figure5 presents temperature profiles at various times. Those profiles allow the evolution of zi with time to be evaluated.Figure 5
Temperature vertical profiles for experiment A2 (the stratification profile is the one att=0 s).Figure6 displays trajectories reconstructed by FT for experiment A2 (Table 1) inside both the stable and unstable layers as they evolve with time. Tracer positions and corresponding trajectories reconstructed over 250 consecutive frames are overlaid, resulting in trajectories extending over a time interval of 10 s. Small segments characterize particles moving with a small velocity while long segments characterize faster particles; identical trajectory lengths correspond to identical velocities and direct picture-to-picture comparisons can be made.Figure 6
Trajectories reconstructed by FT and vertical and horizontal velocity component variance profiles (diamond foru′2¯, triangle for w′2¯) compared with the mixing layer height derived quantitatively from temperature measurement and qualitatively from trajectory observations (experiment A2).


### 3.2. Convective Boundary Layer Growth

The stable boundary layer established during the night is characterized by a positive vertical temperature gradientΓ>0. Figure 4(a) presents the schematic representation of the temperature profile at different times during the convective boundary layer evolution. For low temperature gradients, particles reaching height ze, where buoyancy forces are null, with a nonzero kinetic energy further increase in their height (overshoot). The less the vertical temperature gradient is, the greater the overshoot is observed. Large temperature gradients were employed in the laboratory experiments presented here; the overshoot is then harder to be observed.Further simplifications have been made in order to formalize a bulk model. In Figure4(b), setting zs=0 and zf=zi, a temperature jump occurs at zf=zi (jump model [50]). In Figure 4(c) both the superficial and penetrative layers are neglected (zs=0 and zf=ze=zi). With these assumptions (thermodynamic model [51]), the energy balance equation for a volume of unit base and height zi is:(1)ziρcpdθ=ρcpw′θ′¯sdt-ρcpw′θ′¯z=zidt;being w′θ′¯z=zi=-Aw′θ′¯s, it is(2)dzidt=1+Aw′θ′¯sΓzi.Then(3)zit2=zi2t1+21+Aw′θ′¯st2-t1Γ,where w′θ′¯s is the mean value of w′θ′¯s between t1 and t2.Figure7 presents zi(t) computed for both experiments through:(i)
relation (3), assuming the heat flux is constant;(ii)
temperature: knowing the mean temperature within the mixing layer, the height can be calculated throughzi(t)=(1/Γ)θ-(t)-(θs)t=0,(θs)t=0 being the surface potential temperature at t=0;(iii)
turbulent intensity;(iv)
heat flux derived from the vertical velocity variance.Figure 7
Comparison of CBL growth for experiments A1 and A2 detected with different methods: temperature measurements, velocity data, and relation (3).


### 3.3. Velocity Measurements

Assuming horizontal homogeneity (this hypothesis has been tested and results are reported in [48]), vertical profiles of the horizontal and vertical velocity component variance, σu2 and σw2, respectively, are computed (Figures 8(a) and 8(b)). Results are fairly similar to data published in the literature. Horizontal velocity variances present a uniform profile over the mixing layer and two maxima, one located close to the lower boundary while the other close to the inversion layer. The maxima result from the strong horizontal outflow from the top of thermals as they approach the interface and from roughness in the layer close to the ground. Small variations in the location and value of the upper profile maximum can occur due to the strength of the inversion (i.e., the initial stratification). CBL topped by strong inversions, for instance, produces large horizontal velocity variances because a greater percentage of the upward motion in the thermals is transferred to the horizontal velocity components [31].Figure 8
Time evolution of (a) the horizontal velocity variance profiles; (b) the vertical velocity variance profiles; (c)σw3/z profiles for experiment A2.
(a)
(b)
(c)The vertical velocity variance exhibits a maximum around1/3zi while the inflection point is located approximately at zi. Moreover inside a well-mixed layer where mechanical production can be neglected, the σw3/z versus z profile provides an approximation for the local heat flux, which linearly decreases with height vanishing by definition at z=zi [36], as shown in Figure 8(c). Thus, the extrapolation to zero of the linear regression of σw3/z profiles provides the CBL growth with time.


### 3.4. Lagrangian Integral Scale

The Lagrangian mean velocity, at three different levels, is presented in Figure9 as a function of the nondimensional time, where time is normalized through t*=zi/w*. It is apparent that, although the phenomenon is steady from the Eulerian point of view, this is not true in a Lagrangian framework. The mean velocity varies with time because this statistic is evaluated along the trajectories of particles that move around within the mixing layer, meeting different conditions [52]. Only the asymptotic limit is zero, as long as the velocities are completely uncorrelated. The mean velocity of particles starting close to the ground increases for short times because they cross layers with higher velocity, and an opposite trend is observed for the particles starting close to the capping inversion.Figure 9
Mean vertical velocity component as a function of time in the Lagrangian reference frame at three nondimensional heights.The nonsymmetry of the Lagrangian autocorrelation coefficient of the vertical velocity component, plotted in Figure10 at three different heights, confirms the Lagrangian nonsteadiness of the phenomenon. Since normalization of the correlation coefficients is done by means of the Eulerian variance at the starting level of the trajectories, they are not defined as less than one as in the case of classical Eulerian coefficients; nevertheless in this way the original shape of the function is preserved. The different behavior of the hot updrafts and cold downdrafts can be shown by conditional sampling. The ensemble of trajectories has been divided into two sets: the first one includes trajectories that initially move upwards and the second one trajectories that initially move downwards. Statistics have been computed separately for the two sets (Figure 10). The correlation coefficient of the vertical velocity is not maximum for zero time lag: particles rising from the lower part of the layer (z/zi=0.25) exhibit a maximum correlation for a slightly positive time lag since they are slow in the proximity of the surface and accelerate through the whole depth of the layer. On the contrary, descending particles found in that region decelerate as they get closer to the lower surface and exhibit the maximum correlation at a negative time lag. This behavior is not apparent at midheight (z/zi=0.55) where no sudden changes in particle velocity are imposed by a nearby boundary. In the upper part of the layer (z/zi=0.75) a symmetrical behavior is observed. Due to the presence of the capping inversion, ascending particles decelerate and their correlation is maximum for slightly negative time lags, whereas descending ones accelerate and have a maximum correlation with a positive time lag.Figure 10
Lagrangian autocorrelation coefficient of the vertical velocity component,ρzz, at three nondimensional heights; triangles indicate the updrafts, squares the downdrafts, and solid lines nonconditional sampling analysis.
(a)
(b)
(c)According to this behavior, the1/e Lagrangian time scale, TLa, of the vertical velocity component, computed for positive time lags (Figure 11), decreases with height for the updrafts since this scale is related to the time taken by the particles to reach the capping inversion. The behavior of the descending particles time scale is different since in this case the scale is related to the time that particles take to reach the lower surface. The downdraft time scales are generally larger than those of the updrafts because the former are characterized by a lower velocity, according to the probability density distribution of the vertical velocity component shown in Figure 12.Figure 11
1
/
e Lagrangian time scale of the updraughts (full curve) and downdraughts (broken curve) as a function of the nondimensional height.Figure 12
Probability density distribution of the nondimensional vertical velocity at three nondimensional heights.


### 3.5. Horizontal Spatial Structure

The spatial correlation coefficient(4)ρij(x,r,t)=ui(x,t)uj(x+r,t)¯σiσjassuming the horizontal homogeneity is a function of height z and of the modulus of r, that is, ρij(z,r,t) [48].Figure13 presents the spatial correlation of the vertical component of the velocity field for three instants of time. Velocity vectors belonging to a layer centered around the convective region centerline are taken into account. For small times, the number of velocity samples belonging to the layer within the mixing region taken into account to compute the correlation is inadequate to gather a statistically acceptable result. For large times, the structure dimension becomes comparable to the test section side and the wall effects are not negligible. This then suggests that the velocity data analysis is meaningful from 150 s to about 1000 s from the beginning of the experiments.Figure 13
Spatial autocorrelation of the vertical component of the velocity field for three times (experiment A2).
(a)
(b)
(c)The oscillating behavior of each line allows the transverse dimension of the plumes within the mixing layer to be evaluated for a given time, to which the correlation refers. The distance between two peaks in each plot is compared with the distance between two domes. The characteristic dimension of those structures increases with time. The correlation goes to zero faster for smaller times, in agreement with the expected smaller dimension of the structure; moreover, the dome horizontal characteristic dimension increases with time. The comparison between the results of the correlation analysis and trajectory patterns shows a good agreement for each time instant presented.


### 3.6. Transilient Matrix

A more detailed description of the mixing can be obtained from the analysis of the transilient matrix [1, 52]. Assuming horizontal homogeneity and dividing the investigation field into layers, this matrix represents the fraction of fluid advected from one layer to another in a given time interval. Its first index shows the destination level of mixing and the second the source level. In Figure 14 the transilient matrix is presented for six time intervals; the rows of the matrix are plotted in reverse order to have the height axis increasing upwards. The elements of the secondary diagonal indicate the fraction of fluid that does not change layer in the given time. For small time intervals, the matrix presents values close to one mainly in this region since particles have not had time to mix. As time goes on, larger and larger zones of the matrix are affected by mixing and have nonzero values; the entire convective layer (z/zi<1) is mixed whereas sources and destinations higher than zi remain concentrated close to the secondary diagonal, showing that the capping inversion does not participate in mixing phenomena and pollutants are confined into the convective layer.Figure 14
Transilient matrices for various nondimensional time lags (experiment A2).


### 3.7. Internal Waves

Figure15 shows trajectories in the stable layer. The interaction of the domes with the stable layer produces an oscillation of tracer particles belonging to the stable fluid (internal waves). The oscillatory movement starts in the fluid volume close to the mixing layer and moves upward with time. The same particles with oscillating trajectories are later entrained inside the mixing fluid of increasing height. Internal waves will then be present over the entire stable layer. The oscillation amplitude increases with time, increasing the energy input from the mixing layer. To test whether the measured velocities are due to internal gravity waves, Figure 16 presents the Lagrangian correlation coefficient of the vertical velocity component. The distance between two consecutive peaks of the correlation coefficient provides the internal wave oscillation period, 2π/ω. The frequency of the internal waves, ω, is in a good agreement with the Brunt-Väissälä frequency, N, in both experiments presented.Figure 15
Trajectories reconstructed by FT within the stable layer (case A2; (a) 200 s; (b) 560 s; (c) 600 s; (d) 800 s).
(a)
(b)
(c)
(d)Figure 16
Correlation coefficient forw for case A2.


### 3.8. Comparison with Field Results and Numerical Simulations

Figure17 shows the horizontal and vertical velocity component variance profiles normalized according to the Deardorff mixed layer similarity [38], by employing the square of the convective velocity, w*, and the mixing layer height, zi, computed through temperature measurements in the mixing layer. After normalization, all profiles collapse on the same curve, losing the dependence on time and particular conditions of the experiment.Figure 17
(a) Horizontal and (b) vertical velocity component variance profiles normalized by the Deardorff mixed layer similarity [38] compared with tank experiments and field measurements. The dashed line indicates Adrian et al. [53] results, the full line relation (5).
(a)
(b)The horizontal velocity component variance in the CBL presents a reasonably constant region through the middle of the mixed layer with clear maxima near 0.1zi and 0.75zi. The observation of σw2 normalized profiles proves that the mixing layer upper limit is located where the variance is about 10% of the maximum value, whereas its peak, about 0.5w*2, is in the middle of the mixing layer. At the top and bottom of the CBL, the horizontal velocity component variance is larger than the vertical velocity component one.In Figures17(a) and 17(b), σu2 and σw2 for cases A1 and A2 are compared with laboratory results ([41] (cases S1 and S2) [23, 32, 45, 54]) as well as field measurements (Minnesota and Ashchurch experiments [55] and Phoenix 78 study [56]). These field data have been selected as representative of values reported in the literature for strongly convective conditions [23].The horizontal velocity component variance is similar to the one observed in most of the laboratory models although the constant region assumes in Hibberd and Sawford [23] a larger value equal to 0.34w*2 and the upper peak is weaker in Willis and Deardorff [41] results and negligible in Kumar and Adrian [45] data. The field data are so scattered that it is difficult to discern a general trend of similar shape to that indicated by laboratory data. Our data also match LES profiles by Nieuwstadt et al. [11] (not reported) in terms of both maxima location and variance value in the constant region.The comparison ofσw2 shows a fairly nice agreement. Our results generally fall in the center of the spread of field data with a peak of 0.5w*2 in the middle of the mixed layer. The magnitude of this maximum agrees with results from large eddy simulations (LES) of the CBL of Nieuwstadt et al. [11] as well as the location at 0.3–0.4zi.The dashed line in Figure17(b) represents the faired curve through Adrian et al. [53] detailed nonpenetrative convection experiment. In the upper part of the mixed layer, the experimental data remain consistently higher than Adrian et al. [53] results because of the presence of the entrainment zone. The solid line in Figure 17(b) represents the empirical curve [23]:(5)σw2w*2=1.2zzi2/31-0.9zzi,which approximates the variance behavior until about 1.1zi. A small nonzero variance in the stable layer extends to about 1.5zi and is associated with penetration of thermals and their influence extended above the mean mixed-layer height.


## 3.1. Experimental Procedure

The laboratory model consists of a convection chamber (0.40 × 0.40 × 0.41 m3) containing an initially stable, density stratified fluid, which is heated from below to trigger destabilization and penetrative convection. The fluid stratification is achieved through the two-tank methods [49]. The fluid initial conditions are velocity equal to zero and linear increase of temperature with height through a slope Γ. After the fluid stratification within the tank is completed, a thermostatically controlled hot water bath (set to a temperature Tc) is connected to the bottom metal plate and the experiment may begin. The framed area is Δx=0.200 m long (x-axis) and Δz=0.15 m high (z-axis).Thermocouples are placed within the test section along a vertical array of 26 probes to measure vertical profiles and on the lower boundary to test horizontal homogeneity in supplying heat. Figure5 presents temperature profiles at various times. Those profiles allow the evolution of zi with time to be evaluated.Figure 5
Temperature vertical profiles for experiment A2 (the stratification profile is the one att=0 s).Figure6 displays trajectories reconstructed by FT for experiment A2 (Table 1) inside both the stable and unstable layers as they evolve with time. Tracer positions and corresponding trajectories reconstructed over 250 consecutive frames are overlaid, resulting in trajectories extending over a time interval of 10 s. Small segments characterize particles moving with a small velocity while long segments characterize faster particles; identical trajectory lengths correspond to identical velocities and direct picture-to-picture comparisons can be made.Figure 6
Trajectories reconstructed by FT and vertical and horizontal velocity component variance profiles (diamond foru′2¯, triangle for w′2¯) compared with the mixing layer height derived quantitatively from temperature measurement and qualitatively from trajectory observations (experiment A2).


## 3.2. Convective Boundary Layer Growth

The stable boundary layer established during the night is characterized by a positive vertical temperature gradientΓ>0. Figure 4(a) presents the schematic representation of the temperature profile at different times during the convective boundary layer evolution. For low temperature gradients, particles reaching height ze, where buoyancy forces are null, with a nonzero kinetic energy further increase in their height (overshoot). The less the vertical temperature gradient is, the greater the overshoot is observed. Large temperature gradients were employed in the laboratory experiments presented here; the overshoot is then harder to be observed.Further simplifications have been made in order to formalize a bulk model. In Figure4(b), setting zs=0 and zf=zi, a temperature jump occurs at zf=zi (jump model [50]). In Figure 4(c) both the superficial and penetrative layers are neglected (zs=0 and zf=ze=zi). With these assumptions (thermodynamic model [51]), the energy balance equation for a volume of unit base and height zi is:(1)ziρcpdθ=ρcpw′θ′¯sdt-ρcpw′θ′¯z=zidt;being w′θ′¯z=zi=-Aw′θ′¯s, it is(2)dzidt=1+Aw′θ′¯sΓzi.Then(3)zit2=zi2t1+21+Aw′θ′¯st2-t1Γ,where w′θ′¯s is the mean value of w′θ′¯s between t1 and t2.Figure7 presents zi(t) computed for both experiments through:(i)
relation (3), assuming the heat flux is constant;(ii)
temperature: knowing the mean temperature within the mixing layer, the height can be calculated throughzi(t)=(1/Γ)θ-(t)-(θs)t=0,(θs)t=0 being the surface potential temperature at t=0;(iii)
turbulent intensity;(iv)
heat flux derived from the vertical velocity variance.Figure 7
Comparison of CBL growth for experiments A1 and A2 detected with different methods: temperature measurements, velocity data, and relation (3).


## 3.3. Velocity Measurements

Assuming horizontal homogeneity (this hypothesis has been tested and results are reported in [48]), vertical profiles of the horizontal and vertical velocity component variance, σu2 and σw2, respectively, are computed (Figures 8(a) and 8(b)). Results are fairly similar to data published in the literature. Horizontal velocity variances present a uniform profile over the mixing layer and two maxima, one located close to the lower boundary while the other close to the inversion layer. The maxima result from the strong horizontal outflow from the top of thermals as they approach the interface and from roughness in the layer close to the ground. Small variations in the location and value of the upper profile maximum can occur due to the strength of the inversion (i.e., the initial stratification). CBL topped by strong inversions, for instance, produces large horizontal velocity variances because a greater percentage of the upward motion in the thermals is transferred to the horizontal velocity components [31].Figure 8
Time evolution of (a) the horizontal velocity variance profiles; (b) the vertical velocity variance profiles; (c)σw3/z profiles for experiment A2.
(a)
(b)
(c)The vertical velocity variance exhibits a maximum around1/3zi while the inflection point is located approximately at zi. Moreover inside a well-mixed layer where mechanical production can be neglected, the σw3/z versus z profile provides an approximation for the local heat flux, which linearly decreases with height vanishing by definition at z=zi [36], as shown in Figure 8(c). Thus, the extrapolation to zero of the linear regression of σw3/z profiles provides the CBL growth with time.


## 3.4. Lagrangian Integral Scale

The Lagrangian mean velocity, at three different levels, is presented in Figure9 as a function of the nondimensional time, where time is normalized through t*=zi/w*. It is apparent that, although the phenomenon is steady from the Eulerian point of view, this is not true in a Lagrangian framework. The mean velocity varies with time because this statistic is evaluated along the trajectories of particles that move around within the mixing layer, meeting different conditions [52]. Only the asymptotic limit is zero, as long as the velocities are completely uncorrelated. The mean velocity of particles starting close to the ground increases for short times because they cross layers with higher velocity, and an opposite trend is observed for the particles starting close to the capping inversion.Figure 9
Mean vertical velocity component as a function of time in the Lagrangian reference frame at three nondimensional heights.The nonsymmetry of the Lagrangian autocorrelation coefficient of the vertical velocity component, plotted in Figure10 at three different heights, confirms the Lagrangian nonsteadiness of the phenomenon. Since normalization of the correlation coefficients is done by means of the Eulerian variance at the starting level of the trajectories, they are not defined as less than one as in the case of classical Eulerian coefficients; nevertheless in this way the original shape of the function is preserved. The different behavior of the hot updrafts and cold downdrafts can be shown by conditional sampling. The ensemble of trajectories has been divided into two sets: the first one includes trajectories that initially move upwards and the second one trajectories that initially move downwards. Statistics have been computed separately for the two sets (Figure 10). The correlation coefficient of the vertical velocity is not maximum for zero time lag: particles rising from the lower part of the layer (z/zi=0.25) exhibit a maximum correlation for a slightly positive time lag since they are slow in the proximity of the surface and accelerate through the whole depth of the layer. On the contrary, descending particles found in that region decelerate as they get closer to the lower surface and exhibit the maximum correlation at a negative time lag. This behavior is not apparent at midheight (z/zi=0.55) where no sudden changes in particle velocity are imposed by a nearby boundary. In the upper part of the layer (z/zi=0.75) a symmetrical behavior is observed. Due to the presence of the capping inversion, ascending particles decelerate and their correlation is maximum for slightly negative time lags, whereas descending ones accelerate and have a maximum correlation with a positive time lag.Figure 10
Lagrangian autocorrelation coefficient of the vertical velocity component,ρzz, at three nondimensional heights; triangles indicate the updrafts, squares the downdrafts, and solid lines nonconditional sampling analysis.
(a)
(b)
(c)According to this behavior, the1/e Lagrangian time scale, TLa, of the vertical velocity component, computed for positive time lags (Figure 11), decreases with height for the updrafts since this scale is related to the time taken by the particles to reach the capping inversion. The behavior of the descending particles time scale is different since in this case the scale is related to the time that particles take to reach the lower surface. The downdraft time scales are generally larger than those of the updrafts because the former are characterized by a lower velocity, according to the probability density distribution of the vertical velocity component shown in Figure 12.Figure 11
1
/
e Lagrangian time scale of the updraughts (full curve) and downdraughts (broken curve) as a function of the nondimensional height.Figure 12
Probability density distribution of the nondimensional vertical velocity at three nondimensional heights.


## 3.5. Horizontal Spatial Structure

The spatial correlation coefficient(4)ρij(x,r,t)=ui(x,t)uj(x+r,t)¯σiσjassuming the horizontal homogeneity is a function of height z and of the modulus of r, that is, ρij(z,r,t) [48].Figure13 presents the spatial correlation of the vertical component of the velocity field for three instants of time. Velocity vectors belonging to a layer centered around the convective region centerline are taken into account. For small times, the number of velocity samples belonging to the layer within the mixing region taken into account to compute the correlation is inadequate to gather a statistically acceptable result. For large times, the structure dimension becomes comparable to the test section side and the wall effects are not negligible. This then suggests that the velocity data analysis is meaningful from 150 s to about 1000 s from the beginning of the experiments.Figure 13
Spatial autocorrelation of the vertical component of the velocity field for three times (experiment A2).
(a)
(b)
(c)The oscillating behavior of each line allows the transverse dimension of the plumes within the mixing layer to be evaluated for a given time, to which the correlation refers. The distance between two peaks in each plot is compared with the distance between two domes. The characteristic dimension of those structures increases with time. The correlation goes to zero faster for smaller times, in agreement with the expected smaller dimension of the structure; moreover, the dome horizontal characteristic dimension increases with time. The comparison between the results of the correlation analysis and trajectory patterns shows a good agreement for each time instant presented.


## 3.6. Transilient Matrix

A more detailed description of the mixing can be obtained from the analysis of the transilient matrix [1, 52]. Assuming horizontal homogeneity and dividing the investigation field into layers, this matrix represents the fraction of fluid advected from one layer to another in a given time interval. Its first index shows the destination level of mixing and the second the source level. In Figure 14 the transilient matrix is presented for six time intervals; the rows of the matrix are plotted in reverse order to have the height axis increasing upwards. The elements of the secondary diagonal indicate the fraction of fluid that does not change layer in the given time. For small time intervals, the matrix presents values close to one mainly in this region since particles have not had time to mix. As time goes on, larger and larger zones of the matrix are affected by mixing and have nonzero values; the entire convective layer (z/zi<1) is mixed whereas sources and destinations higher than zi remain concentrated close to the secondary diagonal, showing that the capping inversion does not participate in mixing phenomena and pollutants are confined into the convective layer.Figure 14
Transilient matrices for various nondimensional time lags (experiment A2).


## 3.7. Internal Waves

Figure15 shows trajectories in the stable layer. The interaction of the domes with the stable layer produces an oscillation of tracer particles belonging to the stable fluid (internal waves). The oscillatory movement starts in the fluid volume close to the mixing layer and moves upward with time. The same particles with oscillating trajectories are later entrained inside the mixing fluid of increasing height. Internal waves will then be present over the entire stable layer. The oscillation amplitude increases with time, increasing the energy input from the mixing layer. To test whether the measured velocities are due to internal gravity waves, Figure 16 presents the Lagrangian correlation coefficient of the vertical velocity component. The distance between two consecutive peaks of the correlation coefficient provides the internal wave oscillation period, 2π/ω. The frequency of the internal waves, ω, is in a good agreement with the Brunt-Väissälä frequency, N, in both experiments presented.Figure 15
Trajectories reconstructed by FT within the stable layer (case A2; (a) 200 s; (b) 560 s; (c) 600 s; (d) 800 s).
(a)
(b)
(c)
(d)Figure 16
Correlation coefficient forw for case A2.


## 3.8. Comparison with Field Results and Numerical Simulations

Figure17 shows the horizontal and vertical velocity component variance profiles normalized according to the Deardorff mixed layer similarity [38], by employing the square of the convective velocity, w*, and the mixing layer height, zi, computed through temperature measurements in the mixing layer. After normalization, all profiles collapse on the same curve, losing the dependence on time and particular conditions of the experiment.Figure 17
(a) Horizontal and (b) vertical velocity component variance profiles normalized by the Deardorff mixed layer similarity [38] compared with tank experiments and field measurements. The dashed line indicates Adrian et al. [53] results, the full line relation (5).
(a)
(b)The horizontal velocity component variance in the CBL presents a reasonably constant region through the middle of the mixed layer with clear maxima near 0.1zi and 0.75zi. The observation of σw2 normalized profiles proves that the mixing layer upper limit is located where the variance is about 10% of the maximum value, whereas its peak, about 0.5w*2, is in the middle of the mixing layer. At the top and bottom of the CBL, the horizontal velocity component variance is larger than the vertical velocity component one.In Figures17(a) and 17(b), σu2 and σw2 for cases A1 and A2 are compared with laboratory results ([41] (cases S1 and S2) [23, 32, 45, 54]) as well as field measurements (Minnesota and Ashchurch experiments [55] and Phoenix 78 study [56]). These field data have been selected as representative of values reported in the literature for strongly convective conditions [23].The horizontal velocity component variance is similar to the one observed in most of the laboratory models although the constant region assumes in Hibberd and Sawford [23] a larger value equal to 0.34w*2 and the upper peak is weaker in Willis and Deardorff [41] results and negligible in Kumar and Adrian [45] data. The field data are so scattered that it is difficult to discern a general trend of similar shape to that indicated by laboratory data. Our data also match LES profiles by Nieuwstadt et al. [11] (not reported) in terms of both maxima location and variance value in the constant region.The comparison ofσw2 shows a fairly nice agreement. Our results generally fall in the center of the spread of field data with a peak of 0.5w*2 in the middle of the mixed layer. The magnitude of this maximum agrees with results from large eddy simulations (LES) of the CBL of Nieuwstadt et al. [11] as well as the location at 0.3–0.4zi.The dashed line in Figure17(b) represents the faired curve through Adrian et al. [53] detailed nonpenetrative convection experiment. In the upper part of the mixed layer, the experimental data remain consistently higher than Adrian et al. [53] results because of the presence of the entrainment zone. The solid line in Figure 17(b) represents the empirical curve [23]:(5)σw2w*2=1.2zzi2/31-0.9zzi,which approximates the variance behavior until about 1.1zi. A small nonzero variance in the stable layer extends to about 1.5zi and is associated with penetration of thermals and their influence extended above the mean mixed-layer height.


## 4. Slope Flows

Thermal circulation along inclined planes driven by the horizontal temperature difference between air adjacent to a mountain slope and the ambient air at the same altitude over the neighboring plane (or over the valley center) has been studied by many investigators, principally using field experiments [57–59]. Few water tank laboratory-scale studies of slope flows have been carried out. Deardorff and Willis [60] present experiments carried out in the same convection tank employed for the buoyancy-driven turbulent boundary layer tilted at a 10° angle so that mean-flow shear would occur. By photographically tracking nearly neutrally buoyant oil droplets, various statistics of turbulence have been obtained. An estimate of the turbulence-energy balance within the layer is also presented in an attempt to assess the extent to which baroclinic shear production may cause the turbulence intensities, normalized by the square of the mixed-layer convective velocity scale to deviate from values to be expected in the absence of mean wind shear. Mitsumoto [61] presents a laboratory experiment to simulate the slope wind circulation using a temperature controlled-water tank of dimensions 2.00 × 0.60 × 0.28 m3 containing a slope of 30°. To visualize the flow, polyethylene powder was mixed into the fluid as a tracer. The streaks of each tracer element on photographs were analyzed to obtain the distribution of flow velocity. Vertical profiles of temperature have been measured by thermistors. The author schematically describes flow patterns for typical upslope and downslope circulations. The upslope flow constitutes a circulation of thickness comparable to the slope height, while the downslope circulation is restricted in a lower region with a second circulation above it. Princevac and Fernando [62] present a study on the mechanism responsible for the generation of upslope flow via a laboratory experiment carried out in a tank with a cross section of 0.60 × 0.60 m2 and a depth of 0.70 m. A glycerin and water mixture was used as the working fluid for a typical averaged fluid depth of 50 cm. Suspended particles and fluorescent dye were used to visualize the flow and measure the velocity field. The paper deals with the generation of upslope flows, in particular, the critical angle that is necessary to generate an upslope flow in contrast to the genesis of classical up- and downdraft dominated intrinsic to flat terrain.A simplified theory of slope flows was proposed by Prandtl [9], obtaining a solution in a closed form of the linearized momentum equation. Some corrections to the Prandtl model have been proposed by Gutman and Melgarejo [63] and Ye et al. [64]. Prandtl theory is unable to describe the differences among katabatic currents [6] and anabatic currents [65].


### 4.1. Experimental Procedure

Experiments were performed in a rectangular tank with a length of 1.700 m, a height of 0.210 m, and a width of 0.600 m that is open at the top and has a horizontal aluminum surface at the bottom. A schematic of the flow field facility is shown in Figure18 as well as the reference systems and locations of temperature and velocity profiles. In particular, profile A is in the middle of the slope, profile B at the slope bottom, and profile C in the flat area 80 mm away from the slope bottom. These sets of experiments simulate the slope flow assuming it drains into or moves out of an open area such as a plain or large valley. The temperature of the tank bottom and upper surfaces was controlled by two heat exchangers allowing an approximately uniform thermal stratification of the needed intensity in the water tank. A nearly linear vertical stratification profile, reached after about 16 hours, simulated a stably stratified potential temperature distribution. A flat aluminum plate simulating the slope was mounted above the bottom surface. The plate inclination angle was α=20°. Heating during the daytime and cooling during the nighttime experienced by the slope surface were simulated by means of a series of Peltier cells attached to the inclined plate. For the simulation of the katabatic flow, a downward (negative) heat flux Qslope was considered while for the anabatic flows the heat flux was upward (positive). During the experiments aimed at investigating the anabatic currents, the test section bottom surface was heated to initiate and sustain convection.Figure 18
Experimental apparatus: (A) polystyrene sheet, (B) free surface heat exchanger, (C) free surface thermostat, (D) valley model, (E) bottom surface heat exchanger, (F) camera, (G) bottom surface thermostat, (H) rectangular-shaped electric heater, and (I) reference systems and locations of temperature arrays velocity profiles (units are in m).The instantaneous, two-dimensional velocity field was measured with a temporal resolution of 1/10 s. Temperatures were measured using 27 T-type thermocouples mounted on three rakes placed at representative locations in the flow (Figure18). The thermocouples were equally spaced every 0.002 m in the vertical direction from 0 to 0.018 m. A further array of 23 thermocouples, covering the entire water depth, was used for measuring the initial temperature distribution. Disturbances in the velocity field produced by the racks of temperature sensors were minimized by positioning the racks in a plane parallel to the one illuminated by the light source, placed at a distance of nearly 0.050 m. Table 2 presents the values of characteristic quantities varied in the experiments.Table 2
Slope flow experimental parameters. Timet=0 corresponds to beginning of image and temperature acquisition; for experiments C1 and C2 heat is provided along the slope from t=20 s; for experiment C3 the UHI heat flux begins at time t=2 s and heat is provided along the slope from t=120 s; in experiments C2 and C3 the bottom is heated from t=2 s from the beginning of the experiments to start and sustain convection. Γ was the vertical component of the temperature gradient; Qs was the surplus of surface heat flux due to the UHI; Qslope was slope heat flux. The duration of the experiments is roughly 600 s.
Experiment
Experiment description
Γ
=
(
∂
θ
a
/
∂
z
)  (K/m)
Q
slope  (kW/m²)
Q
s  (kW/m²)
C1
Katabatic flow
100
1.00
—
C2
Anabatic flow
100
1.00
—
C3
Anabatic flow-UHI interaction
100
1.00
0.60In order to understand the main features of the circulation established by thermal effects in an initially stably stratified environment, hybrid Lagrangian particle tracking was employed. This allows one to obtain the velocity and acceleration of passive tracer particles as the first and second derivatives of a moving spline function that filters the particle trajectory coordinates. The velocity fields were obtained by interpolating the velocity data on a regular grid. The mesh size is 0.0025 m × 0.0025 m superimposed onto a region of 0.023 × 0.007 m², resulting in 92 × 28 knots. These velocity data were averaged over time intervals of 2 seconds (20 frames) when the aim was to highlight turbulent structures or 20 seconds (200 frames) to visualize the average circulation. When the velocity field averaged over the smaller time interval was similar to the one averaged over the larger time interval, the former is presented. Velocity profiles have been computed by interpolating on a line sparse velocity data output by HLPT with a spatial resolution of 1 mm and a time interval of 20 seconds.


### 4.2. Katabatic Currents

The model proposed by Manins and Sawford [6] for the investigation of katabatic currents assumes an indefinite plane with slope α. Set the origin of the reference system at the top point where the katabatic current origins (Figure 19) and use the following notation:(i)
H: distance from the slope surface in the z direction to a point largely unaffected by the katabatic flow;(ii)
θr: the reference potential temperature, that is, the ambient temperature for z=0;(iii)
d=θa-θ;(iv)
Q=ρcpw′θ′¯: heat flux, equal to zero at z=H;(v)
wH: vertical velocity component at H;(vi)
U: mean velocity in the x direction;(vii)
h: current thickness;(viii)
δ: buoyancy deficit;(ix)
B=g[(RH-Rs)-Qs]/ρcpθr: buoyancy flux assumed constant;(x)
RH-Rs: divergence of radiation R over height H;(xi)
cD=-u′w′¯s/U2: near wall roughness coefficient.Figure 19
Schematic of the model proposed by Manins and Sawford [6] for the investigation of katabatic currents.DefineU, h, δ, S1, S2, and S3 with the following relations:(6)Uh=∫0Hudz;Uh2=∫0Hu2dz;Uhδ=∫0Hugdθrdz,S1h2δ=∫0Hgdθrzdz;S2hδ=∫0Hgdθrdz;wHH-S3wHd=∫0Hwdz.Usually, with the exception of few cases of small practical relevance, the coefficients introduced may be considered constant (similarity assumption):(7)S1=0.5;S2=0.9;S3=1.After integration, the momentum, energy and mass balance equations result:(8)∂∂tUh=-∂∂x(U2h)-∂∂x12S1h2δcos⁡⁡α+S2hδsin⁡α-cDU2,∂∂tS2hδ+UhN2sin⁡α-S2Ecos⁡⁡α+∂∂xUhδ=B,∂∂xUh=EU,where E is the entrainment coefficient, defined as the ratio between the vertical velocity component at the upper limit of the katabatic layer (wH) and the average horizontal component (U) [66], related to the Richardson number through empirical relations [67]:(9)E=-wHU=AS1Ri+Mwith:(i)
A, M constants;(ii)
Ri=hδcos⁡⁡α/U2 Richardson number.In conditions of neutral atmosphere, previous equations provide a solution as a power law with the following boundary conditions:h=0, Uh=0, and δh=0 at x=0. The solution is(10)h=34Ex,U3=xRi,δ=43ExRi-1/3.Changing Manins and Sawford [6] model by introducing the turbulent buoyancy flux at the top of the katabatic layer into the energy balance equation, Princevac et al. [68] suggest a different law for h, valid also for large values of x:(11)h2=h02+Ax2tan⁡α,where h0 is h at x=0 and A is a parameter to be empirically determined.


#### 4.2.1. Results for Katabatic Winds

The experiment starts when the power source connected to the Peltier cells is turned on and the upper face of the slope is cooled(t=0). Shortly after, a shallow layer near the bottom surface starts cooling, causing the water adjacent to the inclined surface to flow down the slope, forming a katabatic (or drainage) current. Prior to the beginning of the experiment, the fluid within the test section is stably stratified.Figure20 shows the velocity fields for experiment C1. The time interval for data averaging does not influence the results. For this reason, velocity maps obtained averaging data over 2 seconds are presented. This implies the system was unable to capture the turbulent structures characterizing the interface between the principal current and the compensating flow. The velocity field maps show a downslope flow of water cooling towards the bottom of the slope. Cold water flows in the same direction along the plain. The compensating flow slowly sliding upon the katabatic flow can be observed as well. It is weaker and nearly horizontal.Figure 20
Mean velocity fields after (a) 60 s (average over 2 s) and (b) 480 s (average over 2 s) from the beginning of experiment C1.
(a)
t
=
59
–
61 s
(b)
t
=
479
–
481 sFigure21 compares the velocity component along the slope (profile A) and the horizontal velocity (profiles B and C) for case C1. The velocity increases with time, quickly for earlier times (from 60 s to 120 s from the beginning of the experiment) and more gradually in the remaining time lag. In this case, a steady state is not reached. The maximum value of the velocity component along the slope on its centre (profile A) reaches 0.002 m/s, which is much greater than the horizontal component at the base of the slope (profile B; maximum velocity value equal to 0.0008 m/s) and at the flat area downhill (profile C; maximum velocity value equal to 0.0009 m/s). The flow thickness in the middle of the slope is essentially constant with time at a value of about 0.008 m. The thickness variations are negligible even for the other two profile locations (the thickness in B is around 0.020 m whereas it is about 0.015 m in C) and are greater at the base of the slope than in the valley centre. The compensating flow is visible upon the katabatic wind, and it is characterized by lower velocity values (0.0004 m/s) and higher thickness than the downslope flow.Figure 21
Velocity profiles A, B, and C at 5 time intervals from the beginning of experiment C1.


#### 4.2.2. Comparison with Published Data

Comparison with the Manins and Sawford “extended hydraulic approach” is presented in Figure22 where S1 and S2 values are plotted as a function of Rib. The present experimental effort yields S1 values slightly lower than 1, S2 values slightly larger than 1, to be compared to S1=0.5 and S2=0.9 provided by Manins and Sawford [6]. Also included in Figure 22 are the field data of Princevac et al. [69] and the laboratory data of Shindler et al. [70] and Fernando et al. [71], obtained for three different slope angles in the absence of thermal stratification. Our results are in good agreement with the laboratory tank experiments of Fernando et al. [71], in particular for the profile factor S2. The agreement with the field experiments is also reasonably good, even though the bulk Richardson number for the laboratory cases was nearly 3 orders of magnitude lower. This suggests that the profile parameters do not significantly depend neither on the ambient stratification nor on the slope angle.Figure 22
Comparison with Manins and Sawford [6] theoretical model.


### 4.3. Anabatic Currents

Hunt [65] suggests a model to study the anabatic currents that is based on the subdivision of the currents into several layers. Referring to Figure 23, considering a flat surface of constant slope α interposed between two horizontal layers and with:(i)
L*=(u*/w*)3h: Monin-Obukhov length,(ii)
u*: friction velocity,(iii)
h: total height of the anabatic layer,(iv)
k: von Kármán constant,(v)
τ: shear stress,(vi)
q: kinematic heat flux,(vii)
z0: roughness,(viii)
z0T: thermal roughness.Figure 23
Schematic of the model proposed by Hunt [65] and Hunt et al. [7] for the investigation of anabatic currents.Assuming(i)
the slope length is limited to ensure that no separation of the boundary layer will occur [72];(ii)
the Boussinesq hypothesis is valid and the thermal forcing due to buoyancy isgβΔθ;(iii)
nonlinear terms are negligible;(iv)
small slope angle, sosin⁡α=α,the momentum equation projected along thex axis and the energy balance equation give [7](12)∂u∂t=-gβΔθα+∂τ∂z,(13)∂θ∂t=-∂q∂zwith Δθ=θ-θa(z=h).Different velocity and temperature behaviour characterize the different layers.First Superficial Layer
S
s
(
z
0
<
z
<
L
*
). Shear stress prevails over buoyancy forces. Taking into account the conditions at the wall:atz=z0,  u=0;  τ=τs=u*2;atz=z0T, q=qs;   θ=θsyields:(14)uz=u*kln⁡zz0,Δθ-Δθs=-qsku*ln⁡zz0T,τ=ku*z∂u∂z.Second Superficial Layer
S
c (L*<z<hc). Schumann [73] and Hunt [65] suggest the value assumed by hc:(15)ln⁡hcz0=110·k.Hunt [65] further suggests the behavior of the other quantities:(16)u(z)≅uM-uM-u(L*)zL*-1/3,Δθ(z)=Δθ(L*)zL*-1/3,τ=w*z4/3h-1/3∂u∂zwith uM the uniform velocity within the middle layer:(17)uM=λuu*,Δθs=λθqsu*with λu≅λθ≅(1/k)ln⁡⁡(L*/z0)≅10.The shear stress at the interface between the superficial layer and the middle layerτMS can be obtained by integrating (12) and neglecting the nonstationary terms:(18)τMS-τs=α∫0hcβΔθdz=-αgβΔθshc.Middle Layer
M
(
h
c
<
z
<
h
-
h
I
). Neglecting buoyancy forces, (12) provides the shear stress at the interface between the middle layer and the inversion layer [65]:(19)τMI-τMS≅huMΔt≅110uMw*;τMI and h can then be derived.Inversion Layer
I
(
h
-
h
I
<
z
<
h
). The layer depth can be obtained from [7]:(20)τMI≅αN2hI2.


#### 4.3.1. Results for Anabatic Winds

In the simulation of the anabatic current the flow is driven by the buoyancy force associated with the positive heat flux imposed at the slope surface. Initially, the fluid within the test section is stably stratified. The experiment starts by heating the plain to initiate and maintain convection and turning on the power source connected to the Peltier cells to warm the upper face of the slope causing an upslope current. Figure24 describes the mean velocity fields relative to case C2. The upslope flow is easily visible in Figures 24(b) and 24(c) and less so in Figure 24(a) which describes the phenomenon when the slope heating has just started and the incoming flow is very slow. The flow is almost horizontal in the valley, while it reaches the top of the slope through eddy structures which move upslope (Figure 24(c)). They form from the onset of convective cells driven upslope by the pressure gradient along the slope [62]. This flow feature can be visualized only when the velocity fields are averaged over small time interval, 2 seconds in this case. Analogous to case C1, the compensating flow can be observed but is stronger and less horizontal.Figure 24
Mean velocity fields after (a) 60 s (average over 2 s); (b) 480 s (average over 20 s); (c) 480 s (average over 2 s) from the beginning of experiment C2.
(a)
t
=
59
–
61 s
(b)
t
=
470
–
490 s
(c)
t
=
479
–
481 sFigure25 compares the velocity component along the slope (profile A) and the horizontal velocity component (profiles B and C) for case C2. A quasisteady state occurs about 120 s after the heating starts. After 120 s, velocity and flow thickness variations are small and oscillatory. In the middle of the slope, the tangential velocity rapidly increases until 120 s and becomes roughly constant until 600 s. The maximum value is 0.0016 m/s while the maximum thickness is about 0.010 m. The same behavior can be observed in the other two profiles where the maximum velocity value is smaller (0.0012 m/s at the slope base and 0.0011 m/s in the flat area downhill) and the current thickness is larger (between 0.026 and 0.022 m at the beginning and at the end of the experiment for the slope base and between 0.026 and 0.019 m in the flat area downhill) while it decreases towards the end of the experiment. The anabatic current appears thicker at the base of the slope than at the slope center, and it is slightly thicker in the valley centre. The countercurrent compensating flow is roughly 0.040 m thick and has a 0.0006 m/s maximum velocity value.Figure 25
Velocity profiles A, B, and C at 5 time intervals from the beginning of experiment C2.Figure26 presents the temperature profile A for the anabatic and the katabatic current. In the anabatic case, the temperature increase principally manifests itself close to the slope surface until about 0.015 m from the lower boundary. The maximum temperature variation is roughly +0.8 K. At larger heights, temperature profiles uniformly shift toward warmer temperatures (maximum temperature increase of 0.35 K). For the katabatic case, the effects of the heat exchange between the bottom surface and the fluid are clear until about 0.008 m. The temperature decrease (reaching the maximum value of roughly 2 K) is abrupt during the initial 120 s, whereas it becomes gradual until the end of the experiment.Figure 26
Temperature profiles in the middle of the slope for the katabatic case (left-hand side) and the anabatic case (right-hand side).
(a)
(b)As expected, the structure of the anabatic current differs from its katabatic counterpart as a result of the larger turbulent activity and vertical mixing. This is well reflected both in the depth of the upflow current, that is nearly twice that of the downslope one, and in the shape of theu and θ profiles, approximately uniform within the bulk of the mixed layer. The complete set of slope flow experiments [74] suggests that the thickness of the anabatic current shows a dependence on the ambient vertical temperature gradient, similar to the downslope case.


#### 4.3.2. Comparison with Published Data

Experimental results do not evidence a subdivision of the slope thickness similar to Hunt et al. [7]; in particular the superficial layers are not easily distinguishable. Adversely, the parameter λu for all the experiments was computed and it ranges between 1 and 2 (Figure 27). This value is in good agreement with ~2 found by Hunt et al. [7] in the laboratory and is not far off λu≅4.2 found by Princevac et al. [69] for the upslope flow in atmosphere. Note that our value agrees reasonably well with λu=3.3 obtained by Chan [75] in a laboratory simulation of upslope flows.Figure 27
Comparison with Hunt et al. [7] theoretical model.


## 4.1. Experimental Procedure

Experiments were performed in a rectangular tank with a length of 1.700 m, a height of 0.210 m, and a width of 0.600 m that is open at the top and has a horizontal aluminum surface at the bottom. A schematic of the flow field facility is shown in Figure18 as well as the reference systems and locations of temperature and velocity profiles. In particular, profile A is in the middle of the slope, profile B at the slope bottom, and profile C in the flat area 80 mm away from the slope bottom. These sets of experiments simulate the slope flow assuming it drains into or moves out of an open area such as a plain or large valley. The temperature of the tank bottom and upper surfaces was controlled by two heat exchangers allowing an approximately uniform thermal stratification of the needed intensity in the water tank. A nearly linear vertical stratification profile, reached after about 16 hours, simulated a stably stratified potential temperature distribution. A flat aluminum plate simulating the slope was mounted above the bottom surface. The plate inclination angle was α=20°. Heating during the daytime and cooling during the nighttime experienced by the slope surface were simulated by means of a series of Peltier cells attached to the inclined plate. For the simulation of the katabatic flow, a downward (negative) heat flux Qslope was considered while for the anabatic flows the heat flux was upward (positive). During the experiments aimed at investigating the anabatic currents, the test section bottom surface was heated to initiate and sustain convection.Figure 18
Experimental apparatus: (A) polystyrene sheet, (B) free surface heat exchanger, (C) free surface thermostat, (D) valley model, (E) bottom surface heat exchanger, (F) camera, (G) bottom surface thermostat, (H) rectangular-shaped electric heater, and (I) reference systems and locations of temperature arrays velocity profiles (units are in m).The instantaneous, two-dimensional velocity field was measured with a temporal resolution of 1/10 s. Temperatures were measured using 27 T-type thermocouples mounted on three rakes placed at representative locations in the flow (Figure18). The thermocouples were equally spaced every 0.002 m in the vertical direction from 0 to 0.018 m. A further array of 23 thermocouples, covering the entire water depth, was used for measuring the initial temperature distribution. Disturbances in the velocity field produced by the racks of temperature sensors were minimized by positioning the racks in a plane parallel to the one illuminated by the light source, placed at a distance of nearly 0.050 m. Table 2 presents the values of characteristic quantities varied in the experiments.Table 2
Slope flow experimental parameters. Timet=0 corresponds to beginning of image and temperature acquisition; for experiments C1 and C2 heat is provided along the slope from t=20 s; for experiment C3 the UHI heat flux begins at time t=2 s and heat is provided along the slope from t=120 s; in experiments C2 and C3 the bottom is heated from t=2 s from the beginning of the experiments to start and sustain convection. Γ was the vertical component of the temperature gradient; Qs was the surplus of surface heat flux due to the UHI; Qslope was slope heat flux. The duration of the experiments is roughly 600 s.
Experiment
Experiment description
Γ
=
(
∂
θ
a
/
∂
z
)  (K/m)
Q
slope  (kW/m²)
Q
s  (kW/m²)
C1
Katabatic flow
100
1.00
—
C2
Anabatic flow
100
1.00
—
C3
Anabatic flow-UHI interaction
100
1.00
0.60In order to understand the main features of the circulation established by thermal effects in an initially stably stratified environment, hybrid Lagrangian particle tracking was employed. This allows one to obtain the velocity and acceleration of passive tracer particles as the first and second derivatives of a moving spline function that filters the particle trajectory coordinates. The velocity fields were obtained by interpolating the velocity data on a regular grid. The mesh size is 0.0025 m × 0.0025 m superimposed onto a region of 0.023 × 0.007 m², resulting in 92 × 28 knots. These velocity data were averaged over time intervals of 2 seconds (20 frames) when the aim was to highlight turbulent structures or 20 seconds (200 frames) to visualize the average circulation. When the velocity field averaged over the smaller time interval was similar to the one averaged over the larger time interval, the former is presented. Velocity profiles have been computed by interpolating on a line sparse velocity data output by HLPT with a spatial resolution of 1 mm and a time interval of 20 seconds.


## 4.2. Katabatic Currents

The model proposed by Manins and Sawford [6] for the investigation of katabatic currents assumes an indefinite plane with slope α. Set the origin of the reference system at the top point where the katabatic current origins (Figure 19) and use the following notation:(i)
H: distance from the slope surface in the z direction to a point largely unaffected by the katabatic flow;(ii)
θr: the reference potential temperature, that is, the ambient temperature for z=0;(iii)
d=θa-θ;(iv)
Q=ρcpw′θ′¯: heat flux, equal to zero at z=H;(v)
wH: vertical velocity component at H;(vi)
U: mean velocity in the x direction;(vii)
h: current thickness;(viii)
δ: buoyancy deficit;(ix)
B=g[(RH-Rs)-Qs]/ρcpθr: buoyancy flux assumed constant;(x)
RH-Rs: divergence of radiation R over height H;(xi)
cD=-u′w′¯s/U2: near wall roughness coefficient.Figure 19
Schematic of the model proposed by Manins and Sawford [6] for the investigation of katabatic currents.DefineU, h, δ, S1, S2, and S3 with the following relations:(6)Uh=∫0Hudz;Uh2=∫0Hu2dz;Uhδ=∫0Hugdθrdz,S1h2δ=∫0Hgdθrzdz;S2hδ=∫0Hgdθrdz;wHH-S3wHd=∫0Hwdz.Usually, with the exception of few cases of small practical relevance, the coefficients introduced may be considered constant (similarity assumption):(7)S1=0.5;S2=0.9;S3=1.After integration, the momentum, energy and mass balance equations result:(8)∂∂tUh=-∂∂x(U2h)-∂∂x12S1h2δcos⁡⁡α+S2hδsin⁡α-cDU2,∂∂tS2hδ+UhN2sin⁡α-S2Ecos⁡⁡α+∂∂xUhδ=B,∂∂xUh=EU,where E is the entrainment coefficient, defined as the ratio between the vertical velocity component at the upper limit of the katabatic layer (wH) and the average horizontal component (U) [66], related to the Richardson number through empirical relations [67]:(9)E=-wHU=AS1Ri+Mwith:(i)
A, M constants;(ii)
Ri=hδcos⁡⁡α/U2 Richardson number.In conditions of neutral atmosphere, previous equations provide a solution as a power law with the following boundary conditions:h=0, Uh=0, and δh=0 at x=0. The solution is(10)h=34Ex,U3=xRi,δ=43ExRi-1/3.Changing Manins and Sawford [6] model by introducing the turbulent buoyancy flux at the top of the katabatic layer into the energy balance equation, Princevac et al. [68] suggest a different law for h, valid also for large values of x:(11)h2=h02+Ax2tan⁡α,where h0 is h at x=0 and A is a parameter to be empirically determined.


### 4.2.1. Results for Katabatic Winds

The experiment starts when the power source connected to the Peltier cells is turned on and the upper face of the slope is cooled(t=0). Shortly after, a shallow layer near the bottom surface starts cooling, causing the water adjacent to the inclined surface to flow down the slope, forming a katabatic (or drainage) current. Prior to the beginning of the experiment, the fluid within the test section is stably stratified.Figure20 shows the velocity fields for experiment C1. The time interval for data averaging does not influence the results. For this reason, velocity maps obtained averaging data over 2 seconds are presented. This implies the system was unable to capture the turbulent structures characterizing the interface between the principal current and the compensating flow. The velocity field maps show a downslope flow of water cooling towards the bottom of the slope. Cold water flows in the same direction along the plain. The compensating flow slowly sliding upon the katabatic flow can be observed as well. It is weaker and nearly horizontal.Figure 20
Mean velocity fields after (a) 60 s (average over 2 s) and (b) 480 s (average over 2 s) from the beginning of experiment C1.
(a)
t
=
59
–
61 s
(b)
t
=
479
–
481 sFigure21 compares the velocity component along the slope (profile A) and the horizontal velocity (profiles B and C) for case C1. The velocity increases with time, quickly for earlier times (from 60 s to 120 s from the beginning of the experiment) and more gradually in the remaining time lag. In this case, a steady state is not reached. The maximum value of the velocity component along the slope on its centre (profile A) reaches 0.002 m/s, which is much greater than the horizontal component at the base of the slope (profile B; maximum velocity value equal to 0.0008 m/s) and at the flat area downhill (profile C; maximum velocity value equal to 0.0009 m/s). The flow thickness in the middle of the slope is essentially constant with time at a value of about 0.008 m. The thickness variations are negligible even for the other two profile locations (the thickness in B is around 0.020 m whereas it is about 0.015 m in C) and are greater at the base of the slope than in the valley centre. The compensating flow is visible upon the katabatic wind, and it is characterized by lower velocity values (0.0004 m/s) and higher thickness than the downslope flow.Figure 21
Velocity profiles A, B, and C at 5 time intervals from the beginning of experiment C1.


### 4.2.2. Comparison with Published Data

Comparison with the Manins and Sawford “extended hydraulic approach” is presented in Figure22 where S1 and S2 values are plotted as a function of Rib. The present experimental effort yields S1 values slightly lower than 1, S2 values slightly larger than 1, to be compared to S1=0.5 and S2=0.9 provided by Manins and Sawford [6]. Also included in Figure 22 are the field data of Princevac et al. [69] and the laboratory data of Shindler et al. [70] and Fernando et al. [71], obtained for three different slope angles in the absence of thermal stratification. Our results are in good agreement with the laboratory tank experiments of Fernando et al. [71], in particular for the profile factor S2. The agreement with the field experiments is also reasonably good, even though the bulk Richardson number for the laboratory cases was nearly 3 orders of magnitude lower. This suggests that the profile parameters do not significantly depend neither on the ambient stratification nor on the slope angle.Figure 22
Comparison with Manins and Sawford [6] theoretical model.


## 4.2.1. Results for Katabatic Winds

The experiment starts when the power source connected to the Peltier cells is turned on and the upper face of the slope is cooled(t=0). Shortly after, a shallow layer near the bottom surface starts cooling, causing the water adjacent to the inclined surface to flow down the slope, forming a katabatic (or drainage) current. Prior to the beginning of the experiment, the fluid within the test section is stably stratified.Figure20 shows the velocity fields for experiment C1. The time interval for data averaging does not influence the results. For this reason, velocity maps obtained averaging data over 2 seconds are presented. This implies the system was unable to capture the turbulent structures characterizing the interface between the principal current and the compensating flow. The velocity field maps show a downslope flow of water cooling towards the bottom of the slope. Cold water flows in the same direction along the plain. The compensating flow slowly sliding upon the katabatic flow can be observed as well. It is weaker and nearly horizontal.Figure 20
Mean velocity fields after (a) 60 s (average over 2 s) and (b) 480 s (average over 2 s) from the beginning of experiment C1.
(a)
t
=
59
–
61 s
(b)
t
=
479
–
481 sFigure21 compares the velocity component along the slope (profile A) and the horizontal velocity (profiles B and C) for case C1. The velocity increases with time, quickly for earlier times (from 60 s to 120 s from the beginning of the experiment) and more gradually in the remaining time lag. In this case, a steady state is not reached. The maximum value of the velocity component along the slope on its centre (profile A) reaches 0.002 m/s, which is much greater than the horizontal component at the base of the slope (profile B; maximum velocity value equal to 0.0008 m/s) and at the flat area downhill (profile C; maximum velocity value equal to 0.0009 m/s). The flow thickness in the middle of the slope is essentially constant with time at a value of about 0.008 m. The thickness variations are negligible even for the other two profile locations (the thickness in B is around 0.020 m whereas it is about 0.015 m in C) and are greater at the base of the slope than in the valley centre. The compensating flow is visible upon the katabatic wind, and it is characterized by lower velocity values (0.0004 m/s) and higher thickness than the downslope flow.Figure 21
Velocity profiles A, B, and C at 5 time intervals from the beginning of experiment C1.


## 4.2.2. Comparison with Published Data

Comparison with the Manins and Sawford “extended hydraulic approach” is presented in Figure22 where S1 and S2 values are plotted as a function of Rib. The present experimental effort yields S1 values slightly lower than 1, S2 values slightly larger than 1, to be compared to S1=0.5 and S2=0.9 provided by Manins and Sawford [6]. Also included in Figure 22 are the field data of Princevac et al. [69] and the laboratory data of Shindler et al. [70] and Fernando et al. [71], obtained for three different slope angles in the absence of thermal stratification. Our results are in good agreement with the laboratory tank experiments of Fernando et al. [71], in particular for the profile factor S2. The agreement with the field experiments is also reasonably good, even though the bulk Richardson number for the laboratory cases was nearly 3 orders of magnitude lower. This suggests that the profile parameters do not significantly depend neither on the ambient stratification nor on the slope angle.Figure 22
Comparison with Manins and Sawford [6] theoretical model.


## 4.3. Anabatic Currents

Hunt [65] suggests a model to study the anabatic currents that is based on the subdivision of the currents into several layers. Referring to Figure 23, considering a flat surface of constant slope α interposed between two horizontal layers and with:(i)
L*=(u*/w*)3h: Monin-Obukhov length,(ii)
u*: friction velocity,(iii)
h: total height of the anabatic layer,(iv)
k: von Kármán constant,(v)
τ: shear stress,(vi)
q: kinematic heat flux,(vii)
z0: roughness,(viii)
z0T: thermal roughness.Figure 23
Schematic of the model proposed by Hunt [65] and Hunt et al. [7] for the investigation of anabatic currents.Assuming(i)
the slope length is limited to ensure that no separation of the boundary layer will occur [72];(ii)
the Boussinesq hypothesis is valid and the thermal forcing due to buoyancy isgβΔθ;(iii)
nonlinear terms are negligible;(iv)
small slope angle, sosin⁡α=α,the momentum equation projected along thex axis and the energy balance equation give [7](12)∂u∂t=-gβΔθα+∂τ∂z,(13)∂θ∂t=-∂q∂zwith Δθ=θ-θa(z=h).Different velocity and temperature behaviour characterize the different layers.First Superficial Layer
S
s
(
z
0
<
z
<
L
*
). Shear stress prevails over buoyancy forces. Taking into account the conditions at the wall:atz=z0,  u=0;  τ=τs=u*2;atz=z0T, q=qs;   θ=θsyields:(14)uz=u*kln⁡zz0,Δθ-Δθs=-qsku*ln⁡zz0T,τ=ku*z∂u∂z.Second Superficial Layer
S
c (L*<z<hc). Schumann [73] and Hunt [65] suggest the value assumed by hc:(15)ln⁡hcz0=110·k.Hunt [65] further suggests the behavior of the other quantities:(16)u(z)≅uM-uM-u(L*)zL*-1/3,Δθ(z)=Δθ(L*)zL*-1/3,τ=w*z4/3h-1/3∂u∂zwith uM the uniform velocity within the middle layer:(17)uM=λuu*,Δθs=λθqsu*with λu≅λθ≅(1/k)ln⁡⁡(L*/z0)≅10.The shear stress at the interface between the superficial layer and the middle layerτMS can be obtained by integrating (12) and neglecting the nonstationary terms:(18)τMS-τs=α∫0hcβΔθdz=-αgβΔθshc.Middle Layer
M
(
h
c
<
z
<
h
-
h
I
). Neglecting buoyancy forces, (12) provides the shear stress at the interface between the middle layer and the inversion layer [65]:(19)τMI-τMS≅huMΔt≅110uMw*;τMI and h can then be derived.Inversion Layer
I
(
h
-
h
I
<
z
<
h
). The layer depth can be obtained from [7]:(20)τMI≅αN2hI2.


### 4.3.1. Results for Anabatic Winds

In the simulation of the anabatic current the flow is driven by the buoyancy force associated with the positive heat flux imposed at the slope surface. Initially, the fluid within the test section is stably stratified. The experiment starts by heating the plain to initiate and maintain convection and turning on the power source connected to the Peltier cells to warm the upper face of the slope causing an upslope current. Figure24 describes the mean velocity fields relative to case C2. The upslope flow is easily visible in Figures 24(b) and 24(c) and less so in Figure 24(a) which describes the phenomenon when the slope heating has just started and the incoming flow is very slow. The flow is almost horizontal in the valley, while it reaches the top of the slope through eddy structures which move upslope (Figure 24(c)). They form from the onset of convective cells driven upslope by the pressure gradient along the slope [62]. This flow feature can be visualized only when the velocity fields are averaged over small time interval, 2 seconds in this case. Analogous to case C1, the compensating flow can be observed but is stronger and less horizontal.Figure 24
Mean velocity fields after (a) 60 s (average over 2 s); (b) 480 s (average over 20 s); (c) 480 s (average over 2 s) from the beginning of experiment C2.
(a)
t
=
59
–
61 s
(b)
t
=
470
–
490 s
(c)
t
=
479
–
481 sFigure25 compares the velocity component along the slope (profile A) and the horizontal velocity component (profiles B and C) for case C2. A quasisteady state occurs about 120 s after the heating starts. After 120 s, velocity and flow thickness variations are small and oscillatory. In the middle of the slope, the tangential velocity rapidly increases until 120 s and becomes roughly constant until 600 s. The maximum value is 0.0016 m/s while the maximum thickness is about 0.010 m. The same behavior can be observed in the other two profiles where the maximum velocity value is smaller (0.0012 m/s at the slope base and 0.0011 m/s in the flat area downhill) and the current thickness is larger (between 0.026 and 0.022 m at the beginning and at the end of the experiment for the slope base and between 0.026 and 0.019 m in the flat area downhill) while it decreases towards the end of the experiment. The anabatic current appears thicker at the base of the slope than at the slope center, and it is slightly thicker in the valley centre. The countercurrent compensating flow is roughly 0.040 m thick and has a 0.0006 m/s maximum velocity value.Figure 25
Velocity profiles A, B, and C at 5 time intervals from the beginning of experiment C2.Figure26 presents the temperature profile A for the anabatic and the katabatic current. In the anabatic case, the temperature increase principally manifests itself close to the slope surface until about 0.015 m from the lower boundary. The maximum temperature variation is roughly +0.8 K. At larger heights, temperature profiles uniformly shift toward warmer temperatures (maximum temperature increase of 0.35 K). For the katabatic case, the effects of the heat exchange between the bottom surface and the fluid are clear until about 0.008 m. The temperature decrease (reaching the maximum value of roughly 2 K) is abrupt during the initial 120 s, whereas it becomes gradual until the end of the experiment.Figure 26
Temperature profiles in the middle of the slope for the katabatic case (left-hand side) and the anabatic case (right-hand side).
(a)
(b)As expected, the structure of the anabatic current differs from its katabatic counterpart as a result of the larger turbulent activity and vertical mixing. This is well reflected both in the depth of the upflow current, that is nearly twice that of the downslope one, and in the shape of theu and θ profiles, approximately uniform within the bulk of the mixed layer. The complete set of slope flow experiments [74] suggests that the thickness of the anabatic current shows a dependence on the ambient vertical temperature gradient, similar to the downslope case.


### 4.3.2. Comparison with Published Data

Experimental results do not evidence a subdivision of the slope thickness similar to Hunt et al. [7]; in particular the superficial layers are not easily distinguishable. Adversely, the parameter λu for all the experiments was computed and it ranges between 1 and 2 (Figure 27). This value is in good agreement with ~2 found by Hunt et al. [7] in the laboratory and is not far off λu≅4.2 found by Princevac et al. [69] for the upslope flow in atmosphere. Note that our value agrees reasonably well with λu=3.3 obtained by Chan [75] in a laboratory simulation of upslope flows.Figure 27
Comparison with Hunt et al. [7] theoretical model.


## 4.3.1. Results for Anabatic Winds

In the simulation of the anabatic current the flow is driven by the buoyancy force associated with the positive heat flux imposed at the slope surface. Initially, the fluid within the test section is stably stratified. The experiment starts by heating the plain to initiate and maintain convection and turning on the power source connected to the Peltier cells to warm the upper face of the slope causing an upslope current. Figure24 describes the mean velocity fields relative to case C2. The upslope flow is easily visible in Figures 24(b) and 24(c) and less so in Figure 24(a) which describes the phenomenon when the slope heating has just started and the incoming flow is very slow. The flow is almost horizontal in the valley, while it reaches the top of the slope through eddy structures which move upslope (Figure 24(c)). They form from the onset of convective cells driven upslope by the pressure gradient along the slope [62]. This flow feature can be visualized only when the velocity fields are averaged over small time interval, 2 seconds in this case. Analogous to case C1, the compensating flow can be observed but is stronger and less horizontal.Figure 24
Mean velocity fields after (a) 60 s (average over 2 s); (b) 480 s (average over 20 s); (c) 480 s (average over 2 s) from the beginning of experiment C2.
(a)
t
=
59
–
61 s
(b)
t
=
470
–
490 s
(c)
t
=
479
–
481 sFigure25 compares the velocity component along the slope (profile A) and the horizontal velocity component (profiles B and C) for case C2. A quasisteady state occurs about 120 s after the heating starts. After 120 s, velocity and flow thickness variations are small and oscillatory. In the middle of the slope, the tangential velocity rapidly increases until 120 s and becomes roughly constant until 600 s. The maximum value is 0.0016 m/s while the maximum thickness is about 0.010 m. The same behavior can be observed in the other two profiles where the maximum velocity value is smaller (0.0012 m/s at the slope base and 0.0011 m/s in the flat area downhill) and the current thickness is larger (between 0.026 and 0.022 m at the beginning and at the end of the experiment for the slope base and between 0.026 and 0.019 m in the flat area downhill) while it decreases towards the end of the experiment. The anabatic current appears thicker at the base of the slope than at the slope center, and it is slightly thicker in the valley centre. The countercurrent compensating flow is roughly 0.040 m thick and has a 0.0006 m/s maximum velocity value.Figure 25
Velocity profiles A, B, and C at 5 time intervals from the beginning of experiment C2.Figure26 presents the temperature profile A for the anabatic and the katabatic current. In the anabatic case, the temperature increase principally manifests itself close to the slope surface until about 0.015 m from the lower boundary. The maximum temperature variation is roughly +0.8 K. At larger heights, temperature profiles uniformly shift toward warmer temperatures (maximum temperature increase of 0.35 K). For the katabatic case, the effects of the heat exchange between the bottom surface and the fluid are clear until about 0.008 m. The temperature decrease (reaching the maximum value of roughly 2 K) is abrupt during the initial 120 s, whereas it becomes gradual until the end of the experiment.Figure 26
Temperature profiles in the middle of the slope for the katabatic case (left-hand side) and the anabatic case (right-hand side).
(a)
(b)As expected, the structure of the anabatic current differs from its katabatic counterpart as a result of the larger turbulent activity and vertical mixing. This is well reflected both in the depth of the upflow current, that is nearly twice that of the downslope one, and in the shape of theu and θ profiles, approximately uniform within the bulk of the mixed layer. The complete set of slope flow experiments [74] suggests that the thickness of the anabatic current shows a dependence on the ambient vertical temperature gradient, similar to the downslope case.


## 4.3.2. Comparison with Published Data

Experimental results do not evidence a subdivision of the slope thickness similar to Hunt et al. [7]; in particular the superficial layers are not easily distinguishable. Adversely, the parameter λu for all the experiments was computed and it ranges between 1 and 2 (Figure 27). This value is in good agreement with ~2 found by Hunt et al. [7] in the laboratory and is not far off λu≅4.2 found by Princevac et al. [69] for the upslope flow in atmosphere. Note that our value agrees reasonably well with λu=3.3 obtained by Chan [75] in a laboratory simulation of upslope flows.Figure 27
Comparison with Hunt et al. [7] theoretical model.


## 5. Urban Heat Islands

A spatial heat source, such as the one associated with an urban heat island (UHI), generates a thermal plume and a related circulation determined by the differences in temperature and consequently density, with respect to the nearby fluid. The growing plume entrains the surrounding fluid.Field experiments on urban boundary layer have been conducted in many cities as Montreal [76], New York City [77], Cincinnati [78], Columbus [79], St. Louis [80–83], Sapporo [84], and Marseille [85]. Logistic difficulties of field experiments have prompted a number of physical simulations in the laboratory [24, 86, 87]. In particular, we will refer to Lu et al. [24, 86] experiments since they present low-aspect-ratio plume (aspect ratio defined as the plume height to diameter ratio) and formulate a model of the plume with appropriate scaling and similarity parameters so that results can be applied to the full-scale heat island. The quantities measured by Lu et al. [86] in their extensive and systematic water-tank study, conducted in the convection tank of Willis and Deardorff [41], were the mean-temperature field, mixing height and heat-island intensity as functions of surface heating rates, heat-island sizes, and ambient temperature gradients.Lu et al. [24] have proposed a simplified model (bulk model) for a heat island with cylindrical symmetry. The flow is assumed to be turbulent and Reynolds number independent. The following quantities are defined:(i)
D: heat island diameter,(ii)
Qs: mean heat flux within the island at the ground surface,(iii)
zi: plume height where the temperature difference between the island center and the background (far from the island where its effect is negligible) is maximum,(iv)
ze: height where the temperature in the plume center coincides with the ambient temperature,(v)
zr: height at the island boundary where the radial velocity is zero,(vi)
wD=gβD(Qs/ρ0cp)1/3: radial (horizontal) velocity scale (different than w* because in its definition D appears instead of zi),(vii)
δθu: temperature difference between the plume center and the ambient temperature at the ground,(viii)
U: mean radial velocity in the lower portion (to zr where ur=0) of the heated layer near the outer edge of the circular heat island,(ix)
W: mean vertical velocity at zr.Assuming the pressure distribution outside the island is hydrostatic, a pressure deviation must then occur within the plume; this deviation can be obtained by projecting the momentum equation in the vertical direction (∝ stands for proportional): (21)Δp∝ρ0βδθugzi.Assuming zr∝zi and using the mass, momentum, and energy balance equations,(22)UDzi∝WD2;U∝gβδθuzi1/2;δθuW∝Qsρ0cp;therefore(23)U∝wD;W∝wD2ND;zi∝wDN.The ratio between mean vertical and horizontal velocity components is proportional to the Froude number:(24)WU=wDND=Fr.Lu et al. [24] introduce the Reynolds number which corresponds to the square root of the Grashof number when the velocity does not appear in the scaling as an independent parameter:(25)Re2=wD2D2ν2=gβδθuD3ν2=Gr.Ultimately, the nondimensional numbers taken into account are Pr, Gr, and Fr; for a given fluid and thus a given Pr and for high values of Gr, the phenomenon will essentially depend only on Fr.The urban heat island can develop in a stable atmosphere during the nighttime, or in an unstable atmosphere during the daytime (Figure28). One of the fundamental characteristics of the UHI dynamics is the similarity between nighttime N and daytime D UHIs [88]. Circulations are comparable if δθuD=δθuN, where δθuD is the sum of the temperature growth after sunrise and the temperature difference between rural and urban environment (Figure 28). The mentioned similarity is of great interest since it leads to the extension of the Lu et al. [86] bulk convection model, valid for nocturnal UHIs, to daytime UHIs. In that case, U and W would be appropriate velocity scales for daytime UHIs too.Figure 28
Schematic diagram defining (a) nocturnal and (b) diurnal urban heat island intensity.
(a)
(b)


### 5.1. Experimental Procedure

Experiments were performed in the same tank used for the slope flow experiments. The surplus of surface heat flux Q0 between the city and its rural environment is simulated by means of a thin (~0.2 mm), rectangular-shaped electric heater (D=0.05 m ×  L=0.30 m) connected to a suitable power supply. The area under investigation is rectangular, lying in the vertical x-z plane and passing through the center of the electric heater. The framed area is Δx=0.200 m long (x-axis) and Δz=0.048 m high (z-axis). The origin of the reference system is set at the heater center.Many experiments were performed [89]. We will present next the results for a subset of two experiments (Table 3).Table 3
Parameters for the UHI experiments, including the observed values of mixing heightzi and UHI intensity δθu=(θUHI axis-θa)z=0; Qs is the surplus of surface heat flux. The heater dimensions and the buoyancy frequency N=(gβdθa/dz)1/2≅0.43 s−1 are the same for all cases. For case B2 the surface heat flux was 7.2 · 10−5 K m s−1. The duration of the experiments is roughly 800 s.
Experiment number
Experiment description
Q
s (K m s−1)
U (m/s)
z
i (m)
δ
θ
u (K)
Re
Fr
B1
Nocturnal UHI
1.4 · 10−4
0.0026
0.016
1.4
132
0.121
B2
Diurnal UHI
1.4 · 10−4
0.0030
0.028
4.4
151
0.140


### 5.2. Nocturnal UHI

The temperature of the bottom sections is set atθL=θS=297 K, and the temperature of the upper surface is kept at θU=319 K. The resulting (constant with height) ambient temperature gradient is approximately dθa/dz=70 Km−1, corresponding to a Brunt-Väissälä frequency N of approximately 0.42 s−1.Qualitative and quantitative observations start during the UHI circulation growth after the heater is turned on. Figure29(a) shows the trajectories of tracer particles, that is, flow particles, reconstructed by the FT algorithm for the UHI in the nocturnal case. The trajectories, which correspond to a visualized overlapping of the 80 consecutive positions of the tracer particles, extend over a time interval of 8 s. The colors range from blue, which is associated with the first time shown, to red, which is associated with the last time. Both pictures show the flow in a quasisteady state regime, that is, when the circulation is fully established and sidewall effect is still negligible (about 10 minutes after the beginning of the experiment). Lagrangian data are then used to reconstruct time-averaged Eulerian velocity fields through a resampling procedure. A salient feature is the horizontal temporal plume oscillation, an occurrence typically observed when a fluid is heated from below (natural swaying motion). Experimental data shows that the smaller the Qs, the larger the plume oscillation time period. The averaging time interval Δtaver is set to 300 s to filter out the plume oscillation period and to ensure a substantial number of samples are collected to warrant quantitative analysis of turbulence characteristics.Figure 29
Nocturnal UHI: (a) particle trajectories reconstructed by FT; (b) time-averaged Eulerian velocity vectors and streamlines.
(a)
(b)For a givendθa/dz and D, the time needed to reach nearly steady-state flow conditions depends on the horizontal velocity scale U. For all cases considered, the UHI circulation reaches a quasisteady state a few minutes after the heating is started. This condition persists for about 20 min, when the induced circulation arrives at the sidewalls and starts to be influenced by them. Several tests have been conducted to verify the quasisteady-state condition. It is found that when Δtaver is centered at time instants exceeding t=1000 s from the beginning of the experiment, averages and turbulent statistics do not vary appreciably for all cases.Figure30 shows the mean temperature profiles at various distances from the UHI center for the nocturnal case. Except for the superadiabatic surface layer adjacent to the bottom, inside the middle portion of the UHI (x/D=0, 0.2 and 0.4), temperature does not vary appreciably with height. The boundary layer at x/D=1 remains stably stratified with a temperature profile similar to that imposed for the ambient temperature. Temperature profiles suggest that the resulting UHI circulation is dome shaped and characterized by well-mixed conditions within its central region.Figure 30
Mean temperature profiles as a function of the distance from the UHI center for the nocturnal UHI case.
(a)
(b)
(c)
(d)Figure29(b) shows the Eulerian velocity fields and the streamlines for case B1 which provide a topological description of the fluid dynamic fields associated with the UHI circulation. Nearly 1000 velocity samples belonging to the averaging-time interval Δtaver are employed for each of the 90 × 30 grid cells. The shape of the flow patterns for both cases is similar to the trajectories shown in Figure 29(a) which is consistent with the mean temperature profiles previously discussed. All of the images reveal the typical structure of the circulation associated with the UHI, characterized by a convergent flow in the proximity of the ground, an upward flow at the city center, and a divergent flux in the upper layers.


### 5.3. Diurnal UHI

A stable fluid stratification is established above the land and sea sides by settingTL=TS=297 K and TU=319 K. The temperature at the bottom is then rapidly raised to 304 K (the surface heat flux is estimated to be Q0≅300 Wm−2). As a result, a well-mixed layer forms above the rural areas. The electric heater is finally switched on. The UHI then develops within an unstable environment, assuming the features of a daytime UHI. The horizontal velocity scale U has to be calculated by using Qs+Q0 in contrast to Qs, which is used for nocturnal cases.Figure31(a) shows the trajectories of tracer particles for case B2 after the UHI becomes well established. Figure 31(b) displays the averaged velocity field. Its shape is similar to that observed for the corresponding nocturnal UHI (case B1), even though the daytime UHI is wider and stronger relative to its nocturnal counterpart. Comparison among Figures 29 and 31 highlights the differences between the nocturnal and diurnal cases, in terms of both plume height and the presence of convective structures outside the UHI core. The increased UHI depth for the daytime case (zi=0.028 m), compared with the nocturnal one (zi=0.016 m), is related to the presence of the well-mixed layer, which facilitates the vertical development of the thermal plume.Figure 31
Diurnal UHI: (a) particle trajectories reconstructed by FT; (b) time-averaged Eulerian velocity vectors and streamlines.
(a)
(b)Figure32 shows the mean temperature profiles at various distances from the UHI center for the diurnal cases. Both figures illustrate that, except for the superadiabatic surface layer adjacent to the bottom, temperature does not vary appreciably with height. In the diurnal case, temperature does not vary appreciably with height also outside the heater (x/D=1) as well as inside the middle portion of the UHI (x/D=0, 0.2 and 0.4).Figure 32
Mean temperature profiles as a function of the distance from the UHI center for the diurnal UHI case.
(a)
(b)
(c)
(d)


### 5.4. Comparison with Published Data

The ratio of minimum plume width (diameter) to the heat-island diameter, that is, the contraction ratio, falls in the range of 0.2–0.5 reported by other authors [86, 90, 91].Figure33 reports the UHI aspect ratio (zi/D) versus the Froude number for present laboratory results, compared with literature data: numerical results of Catalano et al. [92], Hidalgo et al. [93], Kristóf et al. [94], Kurbatskii [95], Richiardone and Brusasca [96], and Yoshikado [88]; experimental works by Catalano et al. [92], Cenedese and Monti [97], Faust [98], and Lu et al. [24]; field observations of Clarke and McElroy [79]. For Yoshikado [88], values of the surface heat flux and Froude and Reynolds numbers are not reported and have therefore been estimated. Given D=25000 m and dθa/dz=0.007 K m−1, the city radius, and the ambient temperature gradient selected by Yoshikado [88] and assuming H0≅9.6·10-6 K m s−1 as a typical value for nocturnal heat fluxes, the velocity scale, Froude number, and Reynolds number are wD≅3 m s−1, Fr=0.009, and Re=4.8·109. The ensemble of data which appears to collapse on the line of equation zi/D=2.86Fr found by Lu et al. [24] further confirms the validity of the similarity for diurnal and nocturnal cases and providing convincing evidence for the Reynolds-number independence criterion.Figure 33
Ratio of mixing height to diameter(zi/D) as a function of Fr for present laboratory results compared with literature data.To verify the similarity between present numerical and laboratory results and literature data, the normalized vertical profiles of the temperature anomaly and the horizontal and vertical velocity components are compared. Figure34 compares nondimensional temperature profiles in correspondence to the UHI’s centre (x/D=0) with laboratory results by Lu et al. [24] for Re=2920 and Fr=0.089, data from field experiments in Cincinnati, Ohio (Re=1.2·109 and Fr=0.013 [78]), and numerical results by Yoshikado [88] and Catalano et al. [92]. Despite the large differences in Re and Fr, the agreement among the curves is reasonable. This confirms that the shape of the nondimensional mean temperature profiles of low-aspect-ratio UHIs depends on x/D rather than on Reynolds and Froude numbers. Figure 35 compares the profiles of the horizontal velocity component, normalized by U, at x/D=0.5. The agreement among the curves is reasonable; this supports the possibility to extend the Lu et al. [24] theory for the daytime regime. It is remarkable that the maximum intensities for real scale numerical results, including literature data, are larger than those of the laboratory scale data, evidencing a dependence on the fluid type (air, water). This suggests that velocity scale U could be not an appropriate scaling parameter. Further investigation is needed to address this theoretical aspect. Figure 36 reports the profiles for the vertical velocity, normalized by W. Present numerical and laboratory results show a fair agreement with literature data. The maximum vertical velocity is generally attained at 0.5zi.Figure 34
Vertical profiles of the nondimensional temperature anomaly atx/D=0 for the numerical simulations and the laboratory experiments, compared with literature data: numerical results of Catalano et al. [92], Kristóf et al. [94], Kurbatskii [95], and Yoshikado [88]; laboratory experiments of Cenedese and Monti [97] and Lu et al. [24]; field observations of Clarke [78].Figure 35
Normalized horizontal velocity component vertical profiles atx/D=0.5 for the numerical simulations and the laboratory experiments, compared with literature data: numerical results of Catalano et al. [92], Hidalgo et al. [93], Kristóf et al. [94], Kurbatskii [95], Richiardone and Brusasca [96], and Yoshikado [88]; laboratory experiments of Cenedese and Monti [97] and Lu et al. [86].Figure 36
Normalized vertical velocity component vertical profiles atx/D=0 for the numerical simulations and the laboratory experiments, compared with literature data: numerical results of Catalano et al. [92], Hidalgo et al. [93], Kristóf et al. [94], Kurbatskii [95], and Yoshikado [88]; laboratory experiments of Cenedese and Monti [97] and Lu et al. [86].Kristóf et al. [94] present results of a LES simulation in comparison to experimental data presented in Figures 35 and 36. According to qualitative comparisons, the simulation results for the time averaged velocity field are in agreement with the measured vector plot. The width and maximum height of the thermal plume, as well as the direction and the magnitude of velocity vectors, are very similar. The horizontal component along a vertical line at the periphery of the heated plate corresponds closely with the measured values. The measured vertical velocity component along the axis of the heater is slightly overestimated by the simulation, but the tendency and the location of the maximum elevation are well predicted. Measured dimensionless temperature profiles have been remarkably reproduced by the simulation, cross-validating experimental and numerical data.


## 5.1. Experimental Procedure

Experiments were performed in the same tank used for the slope flow experiments. The surplus of surface heat flux Q0 between the city and its rural environment is simulated by means of a thin (~0.2 mm), rectangular-shaped electric heater (D=0.05 m ×  L=0.30 m) connected to a suitable power supply. The area under investigation is rectangular, lying in the vertical x-z plane and passing through the center of the electric heater. The framed area is Δx=0.200 m long (x-axis) and Δz=0.048 m high (z-axis). The origin of the reference system is set at the heater center.Many experiments were performed [89]. We will present next the results for a subset of two experiments (Table 3).Table 3
Parameters for the UHI experiments, including the observed values of mixing heightzi and UHI intensity δθu=(θUHI axis-θa)z=0; Qs is the surplus of surface heat flux. The heater dimensions and the buoyancy frequency N=(gβdθa/dz)1/2≅0.43 s−1 are the same for all cases. For case B2 the surface heat flux was 7.2 · 10−5 K m s−1. The duration of the experiments is roughly 800 s.
Experiment number
Experiment description
Q
s (K m s−1)
U (m/s)
z
i (m)
δ
θ
u (K)
Re
Fr
B1
Nocturnal UHI
1.4 · 10−4
0.0026
0.016
1.4
132
0.121
B2
Diurnal UHI
1.4 · 10−4
0.0030
0.028
4.4
151
0.140


## 5.2. Nocturnal UHI

The temperature of the bottom sections is set atθL=θS=297 K, and the temperature of the upper surface is kept at θU=319 K. The resulting (constant with height) ambient temperature gradient is approximately dθa/dz=70 Km−1, corresponding to a Brunt-Väissälä frequency N of approximately 0.42 s−1.Qualitative and quantitative observations start during the UHI circulation growth after the heater is turned on. Figure29(a) shows the trajectories of tracer particles, that is, flow particles, reconstructed by the FT algorithm for the UHI in the nocturnal case. The trajectories, which correspond to a visualized overlapping of the 80 consecutive positions of the tracer particles, extend over a time interval of 8 s. The colors range from blue, which is associated with the first time shown, to red, which is associated with the last time. Both pictures show the flow in a quasisteady state regime, that is, when the circulation is fully established and sidewall effect is still negligible (about 10 minutes after the beginning of the experiment). Lagrangian data are then used to reconstruct time-averaged Eulerian velocity fields through a resampling procedure. A salient feature is the horizontal temporal plume oscillation, an occurrence typically observed when a fluid is heated from below (natural swaying motion). Experimental data shows that the smaller the Qs, the larger the plume oscillation time period. The averaging time interval Δtaver is set to 300 s to filter out the plume oscillation period and to ensure a substantial number of samples are collected to warrant quantitative analysis of turbulence characteristics.Figure 29
Nocturnal UHI: (a) particle trajectories reconstructed by FT; (b) time-averaged Eulerian velocity vectors and streamlines.
(a)
(b)For a givendθa/dz and D, the time needed to reach nearly steady-state flow conditions depends on the horizontal velocity scale U. For all cases considered, the UHI circulation reaches a quasisteady state a few minutes after the heating is started. This condition persists for about 20 min, when the induced circulation arrives at the sidewalls and starts to be influenced by them. Several tests have been conducted to verify the quasisteady-state condition. It is found that when Δtaver is centered at time instants exceeding t=1000 s from the beginning of the experiment, averages and turbulent statistics do not vary appreciably for all cases.Figure30 shows the mean temperature profiles at various distances from the UHI center for the nocturnal case. Except for the superadiabatic surface layer adjacent to the bottom, inside the middle portion of the UHI (x/D=0, 0.2 and 0.4), temperature does not vary appreciably with height. The boundary layer at x/D=1 remains stably stratified with a temperature profile similar to that imposed for the ambient temperature. Temperature profiles suggest that the resulting UHI circulation is dome shaped and characterized by well-mixed conditions within its central region.Figure 30
Mean temperature profiles as a function of the distance from the UHI center for the nocturnal UHI case.
(a)
(b)
(c)
(d)Figure29(b) shows the Eulerian velocity fields and the streamlines for case B1 which provide a topological description of the fluid dynamic fields associated with the UHI circulation. Nearly 1000 velocity samples belonging to the averaging-time interval Δtaver are employed for each of the 90 × 30 grid cells. The shape of the flow patterns for both cases is similar to the trajectories shown in Figure 29(a) which is consistent with the mean temperature profiles previously discussed. All of the images reveal the typical structure of the circulation associated with the UHI, characterized by a convergent flow in the proximity of the ground, an upward flow at the city center, and a divergent flux in the upper layers.


## 5.3. Diurnal UHI

A stable fluid stratification is established above the land and sea sides by settingTL=TS=297 K and TU=319 K. The temperature at the bottom is then rapidly raised to 304 K (the surface heat flux is estimated to be Q0≅300 Wm−2). As a result, a well-mixed layer forms above the rural areas. The electric heater is finally switched on. The UHI then develops within an unstable environment, assuming the features of a daytime UHI. The horizontal velocity scale U has to be calculated by using Qs+Q0 in contrast to Qs, which is used for nocturnal cases.Figure31(a) shows the trajectories of tracer particles for case B2 after the UHI becomes well established. Figure 31(b) displays the averaged velocity field. Its shape is similar to that observed for the corresponding nocturnal UHI (case B1), even though the daytime UHI is wider and stronger relative to its nocturnal counterpart. Comparison among Figures 29 and 31 highlights the differences between the nocturnal and diurnal cases, in terms of both plume height and the presence of convective structures outside the UHI core. The increased UHI depth for the daytime case (zi=0.028 m), compared with the nocturnal one (zi=0.016 m), is related to the presence of the well-mixed layer, which facilitates the vertical development of the thermal plume.Figure 31
Diurnal UHI: (a) particle trajectories reconstructed by FT; (b) time-averaged Eulerian velocity vectors and streamlines.
(a)
(b)Figure32 shows the mean temperature profiles at various distances from the UHI center for the diurnal cases. Both figures illustrate that, except for the superadiabatic surface layer adjacent to the bottom, temperature does not vary appreciably with height. In the diurnal case, temperature does not vary appreciably with height also outside the heater (x/D=1) as well as inside the middle portion of the UHI (x/D=0, 0.2 and 0.4).Figure 32
Mean temperature profiles as a function of the distance from the UHI center for the diurnal UHI case.
(a)
(b)
(c)
(d)


## 5.4. Comparison with Published Data

The ratio of minimum plume width (diameter) to the heat-island diameter, that is, the contraction ratio, falls in the range of 0.2–0.5 reported by other authors [86, 90, 91].Figure33 reports the UHI aspect ratio (zi/D) versus the Froude number for present laboratory results, compared with literature data: numerical results of Catalano et al. [92], Hidalgo et al. [93], Kristóf et al. [94], Kurbatskii [95], Richiardone and Brusasca [96], and Yoshikado [88]; experimental works by Catalano et al. [92], Cenedese and Monti [97], Faust [98], and Lu et al. [24]; field observations of Clarke and McElroy [79]. For Yoshikado [88], values of the surface heat flux and Froude and Reynolds numbers are not reported and have therefore been estimated. Given D=25000 m and dθa/dz=0.007 K m−1, the city radius, and the ambient temperature gradient selected by Yoshikado [88] and assuming H0≅9.6·10-6 K m s−1 as a typical value for nocturnal heat fluxes, the velocity scale, Froude number, and Reynolds number are wD≅3 m s−1, Fr=0.009, and Re=4.8·109. The ensemble of data which appears to collapse on the line of equation zi/D=2.86Fr found by Lu et al. [24] further confirms the validity of the similarity for diurnal and nocturnal cases and providing convincing evidence for the Reynolds-number independence criterion.Figure 33
Ratio of mixing height to diameter(zi/D) as a function of Fr for present laboratory results compared with literature data.To verify the similarity between present numerical and laboratory results and literature data, the normalized vertical profiles of the temperature anomaly and the horizontal and vertical velocity components are compared. Figure34 compares nondimensional temperature profiles in correspondence to the UHI’s centre (x/D=0) with laboratory results by Lu et al. [24] for Re=2920 and Fr=0.089, data from field experiments in Cincinnati, Ohio (Re=1.2·109 and Fr=0.013 [78]), and numerical results by Yoshikado [88] and Catalano et al. [92]. Despite the large differences in Re and Fr, the agreement among the curves is reasonable. This confirms that the shape of the nondimensional mean temperature profiles of low-aspect-ratio UHIs depends on x/D rather than on Reynolds and Froude numbers. Figure 35 compares the profiles of the horizontal velocity component, normalized by U, at x/D=0.5. The agreement among the curves is reasonable; this supports the possibility to extend the Lu et al. [24] theory for the daytime regime. It is remarkable that the maximum intensities for real scale numerical results, including literature data, are larger than those of the laboratory scale data, evidencing a dependence on the fluid type (air, water). This suggests that velocity scale U could be not an appropriate scaling parameter. Further investigation is needed to address this theoretical aspect. Figure 36 reports the profiles for the vertical velocity, normalized by W. Present numerical and laboratory results show a fair agreement with literature data. The maximum vertical velocity is generally attained at 0.5zi.Figure 34
Vertical profiles of the nondimensional temperature anomaly atx/D=0 for the numerical simulations and the laboratory experiments, compared with literature data: numerical results of Catalano et al. [92], Kristóf et al. [94], Kurbatskii [95], and Yoshikado [88]; laboratory experiments of Cenedese and Monti [97] and Lu et al. [24]; field observations of Clarke [78].Figure 35
Normalized horizontal velocity component vertical profiles atx/D=0.5 for the numerical simulations and the laboratory experiments, compared with literature data: numerical results of Catalano et al. [92], Hidalgo et al. [93], Kristóf et al. [94], Kurbatskii [95], Richiardone and Brusasca [96], and Yoshikado [88]; laboratory experiments of Cenedese and Monti [97] and Lu et al. [86].Figure 36
Normalized vertical velocity component vertical profiles atx/D=0 for the numerical simulations and the laboratory experiments, compared with literature data: numerical results of Catalano et al. [92], Hidalgo et al. [93], Kristóf et al. [94], Kurbatskii [95], and Yoshikado [88]; laboratory experiments of Cenedese and Monti [97] and Lu et al. [86].Kristóf et al. [94] present results of a LES simulation in comparison to experimental data presented in Figures 35 and 36. According to qualitative comparisons, the simulation results for the time averaged velocity field are in agreement with the measured vector plot. The width and maximum height of the thermal plume, as well as the direction and the magnitude of velocity vectors, are very similar. The horizontal component along a vertical line at the periphery of the heated plate corresponds closely with the measured values. The measured vertical velocity component along the axis of the heater is slightly overestimated by the simulation, but the tendency and the location of the maximum elevation are well predicted. Measured dimensionless temperature profiles have been remarkably reproduced by the simulation, cross-validating experimental and numerical data.


## 6. Interaction between Slope Flows and Urban Heat Island

We will present next a case of interaction between two local winds. Figure37 describes the velocity fields relative to case C3. The surplus of surface heat flux between the city and its rural surroundings is simulated by means of a thin, rectangular-shaped electric heater built in a 0.30 × 0.05 m² plastic tape connected to a suitable power supply. The initial 120 s are essentially an urban heat island in a stratified environment only. The slope heating starts at t=120 s. The urban heat island plume can be observed in Figure 37(a), while in Figures 37(b) and 37(c) the upslope flow deforms the plume by bending. The plume is then displaced and destroyed. The heat surplus provided at the UHI location then contributes to the recirculation zone localized about 0.025 m to the left of the central position of the electric heater. The remaining flow is organized into eddies sliding along the slope, as in the anabatic case without urban heat island. The velocity component along the slope increases with time starting from t=120 s to 240 s, that is, 120 s after the heating starts. Also in this case the steady state is reached after about 120 s.Figure 37
Mean velocity fields after (a) 60 s (average over 2 s), (b) 480 s (average over 20 s), and (c) 480 s (average over 2 s) from the beginning of case C3.
(a)
t
=
59
–
61 s
(b)
t
=
470
–
490 s
(c)
t
=
479
–
481 sFigure38 compares the velocity component along the slope (profile A) and the horizontal velocity components (profiles B and C) for case C3. Profile A shows that the anabatic current starts after 120 s from the beginning of the experiment. Velocity quickly increases until 240 s and slows down until 480 s reaching 0.002 m/s. The velocity distribution is similar to case C1, with larger velocity values. Profile B clearly shows the rapid increase in the horizontal velocity after the slope heating starts and the interaction with the urban heat island disturbing the horizontal flow directed from the flat area downhill to the slope. The profile, significantly different than Figure 25, shows the effects of turbulent mixing which tends to homogenize quantities in the vertical direction. After 240 s, the velocity profile has an almost uniform value roughly equal to 0.0007 m/s, for a depth of 0.015 m. Higher velocity linearly decreases until −0.0005 m/s. The velocity is zero at a distance greater than 0.020 m from the lower boundary. Velocity increases quickly until 240 s and moderately until 480 s. The flow thickness does not change substantially in time both at the slope base and center, but its value is much greater at the base, while its velocity is much lower. The effects of the anabatic wind on the urban heat island can be summarized as follows (profile C): the horizontal velocity component is almost zero for the initial 60 s; then (until 180 s) the profile changes assuming the characteristic shape of the velocity distribution expected within a plume, close to its center; after 180 s, when the anabatic current is well established, the velocity increases and the profile shape resembles the one reconstructed along the slope with a comparable maximum velocity value. The anabatic layer thickness is more than 0.010 m. The return current is similar to the one reconstructed for case C1.Figure 38
Velocity profiles A, B, and C at 5 time intervals from the beginning of case C3.


## 7. Conclusions

In laboratory simulations of the atmospheric boundary layer the geometrical and dynamical similarity with the real case cannot be reproduced. Reynolds and Rayleigh numbers in laboratory models are inherently lower than those in the natural ABL. The maximum values of Ra that can be achieved are on the order 108 whereas in real atmosphere Ra is several orders of magnitude larger (Ra=1018). Nevertheless, the system can be accurately modeled when the Rayleigh number is sufficiently high (Ra larger than 106), and then the main details of the flow are independent of the Reynolds and Rayleigh numbers. Therefore, when these dimensionless numbers exceed some critical value, the flow field simulated in a water tank becomes similar to that of the natural system [99].The advantage of laboratory models compared to the field surveys is that all parameters can be controlled (in particular the fluid stratification and the heat flux at the boundary) and measurements can be undertaken in all parts of the field. The advantage of laboratory models compared with numerical simulations is that the latter uses models of turbulence closure based mostly on assumptions not fully verified.Velocity and temperature measurements have been made on temperature-controlled hydrodynamic tanks. The absence of geostrophic wind is assumed. Water is employed instead of air to obtain higher Rayleigh numbers for the same tank dimensions. Water also permits the use of tracer particles with the same weight of the working fluid, allowing a multipoint monitoring of the phenomenon time evolution using not intrusive optical methods. The large amount of data available makes it possible to compute robust statistics useful for the characterization of turbulence.The effects of the forcing that cause the local winds have been analyzed separately. It was pointed out as there is no generalized theory of similarity. For every situation, the fundamental parameters that characterize the phenomenon have to be highlighted through empirical or semiempirical evaluations. The identification of these parameters is often done through bulk model using very rough schematizations of the balance equations.The analysis of the time evolution of the phenomenon was studied with particular care, highlighting the different time scales related to turbulence, waves, and daily changes.If possible, laboratory results have been compared with other published results. Those comparisons are generally positive, confirming what has been most often presented in the literature.These phenomena are usually very complicated and their interpretation sometimes is not yet fully convincing. Even more complex and difficult to analyze is the case of the interaction of several causes originating local winds. The interaction has been highlighted with two special cases when a heat island interferes with the katabatic current or a sea breeze.



---

*Source: 618903-2015-03-22.xml*