Abstract

In many micrometeorological studies with computational fluid dynamics, building-resolving models usually assume a neutral atmosphere. Nevertheless, urban radiative transfers play an important role because of their influence on the energy budget. To take into account atmospheric radiation and the thermal effects of the buildings in simulations of atmospheric flow and pollutant dispersion in urban areas, a three-dimensional (3D) atmospheric radiative scheme has been developed in the atmospheric module of the Code_Saturne 3D computational fluid dynamic model. On the basis of the discrete ordinate method, the radiative model solves the radiative transfer equation in a semitransparent medium for complex geometries. The spatial mesh discretization is the same as the one used for the dynamics. This paper describes ongoing work with the development of this model. The radiative scheme was previously validated with idealized cases. Here, results of the full coupling of the radiative and thermal schemes with the 3D dynamical model are presented and are compared with measurements from the Mock Urban Setting Test (MUST) and with simpler modeling approaches found in the literature. The model is able to globally reproduce the differences in diurnal evolution of the surface temperatures of the different walls and roof. The inhomogeneous wall temperature is only seen when using the 3D dynamical model for the convective scheme.

1. Introduction

Interest in urban climatology has increased in the past decade. It corresponds to the thermal and dynamical airflow response to the urban system solicitations, resulting in radiative transfers and convective exchanges within the urban air and with the building walls (Grimmond and Oke 1999; Arnfield 2003). In the past few years, numerical studies have been conducted to solve the surface energy balance (SEB) in urban canopies, with different degrees of simplification, using either an integrated representation of the urban canopy (Masson 2000) or a three-dimensional approach (Mills 1996; Miguet and Groleau 2002; Kanda et al. 2005; Krayenhoff and Voogt 2007; Gastellu-Etchegorry 2008; Asawa et al. 2008). Those models share the following parameterizations in their design: the schemes possess separate energy budgets for roofs, roads, and walls; radiative interactions between roads and walls are explicitly treated.

The Town Energy Balance (TEB) scheme of Masson (2000) consists of a facet-averaged scheme with one generic roof, one generic wall, and one generic road. The advantage of the integrated resolution is that few individual SEBs need to be resolved and therefore computation time is kept low, with a simple approach to model the inner-canopy wind flow. TEB has been shown to reproduce accurately the SEB from regional to mesoscale and urban scales (Masson et al. 2002; Lemonsu et al. 2004). Mills (1996) developed the Urban Canopy-Layer Climate Model, which has a detailed representation of the canyon with a highly simplified wind parameterization. The “SOLENE” (Miguet and Groleau 2002) and 3D-Computer Aided Design (Asawa et al. 2008) models are based on a realistic description of the canopy structure using a geometrical 3D surface model assigning radiative and thermal properties to each subfacet of the model and a constant transfer coefficient for each class of elements. The model is originally designed for simulating sunshade, natural lighting, and heat transfers for architectural purposes. The Discrete Anisotropic Radiative Transfer (DART) model (Gastellu-Etchegorry et al. 2004) simulates the radiative transfer in the whole optical domain simultaneously in the atmosphere and in the urban and vegetated landscapes, with or without topography. A major feature of DART is that it can simulate images in the plane of the sensor, for different altitudes from the bottom to the top of the atmosphere. The new version of the DART model, DART-EB, (Gastellu-Etchegorry 2008) includes an energy balance component. In the case of urban canopies, turbulent fluxes and conduction are computed with classical boundary layer laws, using the equations of the TEB model (Masson 2000). The Simple Urban Energy Balance Model for Mesoscale Simulations (SUMM; Kanda et al. 2005), which represents the urban canopy with an infinitely extended regular array of uniform buildings, is more adapted for the mesoscale. The Temperatures of Urban Facets in 3D (TUF-3D) model (Krayenhoff and Voogt 2007) uses the radiosity approach based on interpatch view factors to model radiative exchange between the identical square patches that compose the simplified 3D urban geometry. An exponential inner-canopy wind speed profile is employed. TUF-3D has applications in both surface temperature distributions and thermal remote sensing anisotropy at several scales.

Previously described models have all put a strong emphasis on radiative exchanges but not on a detailed flow field. In this work, in addition to the above applications, we are also interested in applying the model to pollutant dispersion in low–wind speed conditions, when the thermal effects have a strong influence on the flow.

To model the airflow in the urban canopy in nonneutral conditions more accurately and to take into account the 3D convective exchanges, we developed a 3D microscale radiative model coupled with a 3D computational fluid dynamics (CFD) code for complex geometries to simulate dynamics and thermodynamics of the urban atmosphere (Milliez 2006). Differing from other radiative models that calculate the view factors to estimate the incoming radiative fluxes on urban surfaces, our model directly solves the 3D radiative transfer equation in the whole fluid domain. This approach allows us to determine the radiation flux not only on the facets of the urban landscape but also in each fluid grid cell between the buildings. The difference could become important in the case of smoke or fog between the buildings. The model was evaluated with idealized cases, using as a first step a constant 3D wind field (Milliez et al. 2006). The purpose of the work presented here is to study the full radiative–dynamical coupling, using an evolving 3D flow field. First we present the model, and then we discuss in detail the results of the full coupling. We further discuss the influence on the surface temperature of the internal building temperature and the wall thermal modeling, comparing the 3D resolution with the approaches used in other models.

2. Equations and model design

As a key parameter, surface temperature Tsfc is determined by the SEB, which governs the energy exchange processes between each urban surface and the atmosphere (Fig. 1). It is given by

 
formula

where Qcond is the conductive heat flux (W m−2) within the building or the ground subsurface that links the surface temperature to the internal-building or the deep-soil temperature, QH is the sensible heat flux (W m−2) and depends on the local wind intensity, S* is the net shortwave radiative flux (W m−2), and L* is the net longwave radiative flux (W m−2). We neglect in this study the other energy fluxes such as the anthropogenic flux and the latent heat flux. In our model the advection fluxes are taken into account by the full resolution of the flow field.

Fig. 1.

Schematic representation of the energy exchanges at an urban surface (cross section of a building) (S*: net shortwave radiative flux; L*: net longwave radiative flux; QH: sensible heat flux; Qcond: conduction heat flux; TW: wall surface temperature; Tint: internal building temperature).

Fig. 1.

Schematic representation of the energy exchanges at an urban surface (cross section of a building) (S*: net shortwave radiative flux; L*: net longwave radiative flux; QH: sensible heat flux; Qcond: conduction heat flux; TW: wall surface temperature; Tint: internal building temperature).

a. CFD model

To solve the dynamics and therefore to resolve the QH term explicitly, simulations are performed with the 3D open-source CFD code known as Code_Saturne (Archambeau et al. 2003), which can handle complex geometry and complex physics. The flow features in built-up areas make the modeling within the urban canopy difficult. Some typical effects that we have to handle are 3D vortices behind the buildings, high wind speed near the edges of the upwind face, wake effects, and modified turbulence.

In this work, we use the atmospheric module of Code_Saturne, described in detail in Milliez and Carissimo (2007), which takes into account the larger-scale meteorological conditions and the thermal stratification of the atmosphere. In our simulations, we use a Reynolds-averaged Navier–Stokes (RANS) approach with a kɛ turbulence closure. The numerical solver is based on a finite-volume approach for collocated variables on an unstructured grid. Time discretization is achieved through a fractional step scheme, with a prediction–correction step.

b. Radiative model

A new atmospheric 3D radiative scheme was developed in Code_Saturne for the urban canopy (Milliez 2006). We have adapted to the atmosphere a radiative heat transfer scheme available for complex geometry in Code_Saturne that solves the radiative transfer equation for a gray nondiffusive semitransparent media:

 
formula

where I(x, D) is the intensity of radiation at the point x and for the propagation direction D, K(x) is the absorption coefficient, L(x, D) is the luminance (W m−3 sr−1), and Ta is the air temperature (K). In a semitransparent media, I (W m−3 sr−1) and K can be considered to be independent of the wavelength and are integrated over the spectrum. The rate of radiation heating Srad (W m−3) is then given by

 
formula

where dΩ is the element of the solid angle (sr) around the direction.

1) Discrete ordinate method (DOM)

To solve the radiative transfer equation, we chose the discrete ordinate method (Fiveland 1984; Truelove 1987; Liu et al. 2000), which is based on the directional propagation of the radiative wave. The spatial discretization uses the same mesh as the CFD model. The angular discretization has two resolutions: 32 or 128 directions.

2) Shortwave and longwave radiation

As is usually done, we separate the atmospheric radiation into shortwave and longwave radiation. The total incoming and outgoing shortwave radiative fluxes for each solid surface are given by

 
formula
 
formula

where S and S are respectively the incoming and outgoing shortwave radiative fluxes (W m−2), SD is the direct solar flux (W m−2), Sf is the solar flux diffused by the atmosphere above our simulation domain (W m−2), Se is the flux diffused by the environment, that is, resulting from the multireflections on the other subfacets (W m−2), and α is the albedo of the surface.

We express the longwave radiation flux for each surface as

 
formula
 
formula

where L and L are respectively incoming and outgoing longwave radiation flux (W m−2), ɛ is emissivity of the surface; σ is the Stefan–Boltzmann constant (5.667 03 × 10−8 W m−2 K−4), Tsfc is the surface temperature (K), and La and Le are the longwave radiation flux from the atmosphere and from the multireflection on the other surface. As the first step of validation, we assume that, at the scale of our simulations, the atmosphere between the buildings is transparent and set the absorption coefficient to 0 for both the longwave and shortwave radiation.

3) Surface temperature model

The force–restore approach (Deardorf 1978) is commonly used to calculate the ground temperature in meteorological models. This approach is considered to be a very useful tool because a prognostic equation for temperature is used to reproduce the response to periodic heating of the soil. This model has been extended to urban surfaces (Johnson et al. 1991; Dupont and Mestayer 2006). This extension nevertheless supposes well-insulated buildings with a nearly constant internal temperature and homogeneous material. In our model, the force–restore method has been available for some time in simple geometries and has been extended to complex geometries. Because of the limiting hypotheses built into the method, however, especially concerning the deep-soil temperature, we have also tested a simple wall thermal model with a given thickness and an average thermal conductivity.

(i) Force–restore model

The force–restore model is based on a two-layer decomposition of a material considered to be homogeneous: the surface-layer temperature Tsfc responding to external forcing and the deeper layer independent of the diurnal variation. It reads

 
formula

where ω is the Earth angular frequency (Hz), μ is the thermal admittance (J m−2 s−0.5 K−1), and Tg/b is either deep-soil or internal building temperature (K).

(ii) Wall thermal model

This model solves the conduction equation to compute the wall temperature. It reads, after expressing each term in Eq. (1), as

 
formula

where λ is the average thermal conductivity of the wall (W K−1 m−1), e is the thickness of the wall (m), Tint is the internal air temperature (K), hf is the heat transfer coefficient (W m−2 K−1) computed from local flow parameters, and Ta is the external air temperature (K).

4) Internal building temperature

In a real building with good insulation, the variation of the internal building temperature is small. In the experiment we simulate (see section 3), however, the buildings are represented by poorly insulated shipping containers. In this case, the variation of the internal temperature is important and has a great influence on the surface temperature. The internal temperature was not measured in the experiment, however and we computed it with one of the following methods.

(i) Constant T

In this case, the internal building temperature is set to a constant and is computed by averaging the diurnal temperatures of all of the building surfaces.

(ii) Evolution equation

A temperature evolution equation, as in Masson et al. (2002), is used to represent the temperature inside the buildings:

 
formula

where and are the computed internal temperatures (K) at the following and previous time step, respectively, Δt is the time step (s), τ is the period (equal to 1 day; (s), and is the average over all of the surface temperatures (K) computed at time step n.

(iii) Evolution equation with interpolated from measurements

We use the previous formula [Eq. (10)] and replace from the calculation with the average of the measured surface temperatures.

Figure 2 compares the northeast wall surface temperature of a shipping container computed with the three different internal temperature models and the measurements from the Mock Urban Setting Test experiment (see section 3). The value of the constant average internal temperature model is 24°C, which is approximately 2 times the initial values in the other two internal temperature models. Before sunrise, this high value of the constant internal temperature induces a too-rapid heating of the northeast wall. At midday, the constant internal temperature is too low to account for the warming of the surface by the interior air heated by the other sunlit surfaces of the container. Computing the internal temperature with an evolution equation model decreases the heating of the northeast wall before sunrise and improves the results at midday. The differences in the results obtained by using the two evolution equation models (with from computation and from measurements) are small, which was expected because Δt/τ is small. Nevertheless, using the observed average temperature slightly improves the results, especially during the night, when the atmospheric radiative fluxes decrease. The same conclusion applies for the other sides of the containers, which are not shown here. We stress that this may be relevant only in the case of metal containers.

Fig. 2.

Comparison of the diurnal evolution of the surface temperature of a northeast-oriented wall computed with three different models for the internal temperature (times signs: surface temperature from the measurements; dashed line: taking a constant internal building temperature; solid line: taking an evolution equation; filled circles: taking the average surface temperature from the measurements and then using the result in the evolution equation).

Fig. 2.

Comparison of the diurnal evolution of the surface temperature of a northeast-oriented wall computed with three different models for the internal temperature (times signs: surface temperature from the measurements; dashed line: taking a constant internal building temperature; solid line: taking an evolution equation; filled circles: taking the average surface temperature from the measurements and then using the result in the evolution equation).

5) Convection model

The thermal energy equation of the flow must be solved, both to determine stratification effects on vertical turbulent transport and to estimate the surface–air thermal gradient that controls convective heat transfer. The sensible heat flux QH is given by

 
formula

Detailed comparisons between different approaches to model the heat transfer coefficient will be discussed in section 5. Our CFD model solves in 3D the RANS equations in the entire fluid domain. In our simulations, we use a rough-wall boundary condition. The hf is computed for each solid subfacet, depending on the local friction velocity:

 
formula

where ρ is flow density (kg m−3), Cp is specific heat (J kg−1 K−1), u* is the friction velocity, κ is the von Kármán constant, σt is the turbulent Prandtl number, d is the distance (m) to the wall, z0 is the roughness length (m), is the thermal roughness length (m), and fm and fh are the Louis (1979) stability functions.

3. Study case: The Mock Urban Setting Test

a. Configuration of the experiment

The Mock Urban Setting Test (MUST; Biltoft 2001; Yee and Biltoft 2004) conducted in the Utah desert is a near-full-scale experiment that consists of measurements in an idealized urban area represented by 120 shipping containers (length L × width W × height H = 12.2 m × 2.42 m × 2.54 m) arranged in a regular array. MUST has already been used to validate dynamics and dispersion models (Brook et al. 2002; Hanna et al. 2002; Camelli et al. 2005; Milliez and Carissimo 2007, 2008). Because temperature data are also provided, we used the MUST field experiment to study in detail the dynamic–radiative coupling. We focused our study on one instrumented container within the array, and therefore the computational domain was reduced to three rows of three containers (Fig. 3).

Fig. 3.

Mesh of the domain and the subdomain for the 0.8 m × 0.5 m × 0.5 m resolution.

Fig. 3.

Mesh of the domain and the subdomain for the 0.8 m × 0.5 m × 0.5 m resolution.

From the MUST experiment, we selected the day of 25 September 2001. Despite a fairly strong wind (Umean = 7 m s−1), we selected this day because it has already been partly simulated for studies on dispersion (Milliez and Carissimo 2007) and moreover because a complete 24-h dataset for the upstream wind and the surface temperature was available (which was not the case for other days). During this day (Fig. 4), the wind velocity varied from Umin = 3 m s−1 to Umax = 11.5 m s−1 and the average diurnal air temperature was about 24°C (measured at 10 m). For our coupling study, the wind speed may be a little high to test strong radiative effects on the airflow, but it emphasizes the convective effects on the surface temperature.

Fig. 4.

Meteorological data at 10 m from the MUST experiment for 25 Sep 2001: (a) wind speed (m s−1), (b) meteorological wind direction (°), and (c) wind (air) temperature (°C).

Fig. 4.

Meteorological data at 10 m from the MUST experiment for 25 Sep 2001: (a) wind speed (m s−1), (b) meteorological wind direction (°), and (c) wind (air) temperature (°C).

b. Boundary conditions and surface parameters

Table 1 gives the parameters used in the simulation. To be consistent with the experiment, the wind inlet boundary conditions are determined from measurements, using a meteorological file that contains the wind velocity, turbulence kinetic energy, dissipation rate, and temperature profiles for every 2 h. The variation of the deep-soil temperature is neglected. The internal building temperature is computed by the evolution equation with from measurements. We take the same value of the roughness length z0 as in Eichhorn and Balczo (2008). The thermal roughness length is considered to be 1/10 of z0 (Garratt 1992). Because some thermal properties were not determined during the experiment, their values were taken from the literature: we took the values of the albedo and emissivity of the wall to be those of corrugated iron (Oke 1987). The thermal admittance was taken assuming an insulating material in the walls, as observed in some pictures. The ground albedo αs that was input to the model was evaluated from the incoming and outgoing solar fluxes measured upstream by pyranometers and depends on the zenith angle.

Table 1.

Parameters employed in the MUST simulations, with Zref being the reference height for forcing data, Uref being the initial wind speed at Zref, and θ being the initial wind angle. Soil z0, wall z0, and roof z0 are respectively the roughness length for momentum of soil, walls, and roofs. Soil , wall , and roof are respectively the roughness length for heat of soil, walls, and roofs; αS, αW, and αR are respectively the shortwave albedo of soil, walls, and roofs; ɛS, ɛw, and ɛR are respectively the longwave emissivity of soil, walls, and roofs; μS, μW, and μR are respectively the admittance of soil, walls, and roofs; Tair is the initial air temperature at Zref; TS and TSint are respectively the initial soil temperature and deep-soil temperature; TTop, TNW, TSE, TNE, and TSW are initial individual surface temperature; TWint and TRint are respectively the building and roof internal temperature.

Parameters employed in the MUST simulations, with Zref being the reference height for forcing data, Uref being the initial wind speed at Zref, and θ being the initial wind angle. Soil z0, wall z0, and roof z0 are respectively the roughness length for momentum of soil, walls, and roofs. Soil , wall , and roof  are respectively the roughness length for heat of soil, walls, and roofs; αS, αW, and αR are respectively the shortwave albedo of soil, walls, and roofs; ɛS, ɛw, and ɛR are respectively the longwave emissivity of soil, walls, and roofs; μS, μW, and μR are respectively the admittance of soil, walls, and roofs; Tair is the initial air temperature at Zref; TS and TSint are respectively the initial soil temperature and deep-soil temperature; TTop, TNW, TSE, TNE, and TSW are initial individual surface temperature; TWint and TRint are respectively the building and roof internal temperature.
Parameters employed in the MUST simulations, with Zref being the reference height for forcing data, Uref being the initial wind speed at Zref, and θ being the initial wind angle. Soil z0, wall z0, and roof z0 are respectively the roughness length for momentum of soil, walls, and roofs. Soil , wall , and roof  are respectively the roughness length for heat of soil, walls, and roofs; αS, αW, and αR are respectively the shortwave albedo of soil, walls, and roofs; ɛS, ɛw, and ɛR are respectively the longwave emissivity of soil, walls, and roofs; μS, μW, and μR are respectively the admittance of soil, walls, and roofs; Tair is the initial air temperature at Zref; TS and TSint are respectively the initial soil temperature and deep-soil temperature; TTop, TNW, TSE, TNE, and TSW are initial individual surface temperature; TWint and TRint are respectively the building and roof internal temperature.

c. Numerical sensitivity

Before starting our simulations, we performed a numerical sensitivity study for the grid size, the number of directions in the discrete ordinate method, and the radiative time step. We focus here on the radiative aspects because a sensitivity study of the dynamical part has already been performed (Milliez and Carissimo 2007, 2008).

Figure 5 shows surface temperature evolution with grids of different resolutions. The fine (~55 000 cells) and the refined grids (~173 000 cells) give similar results, and the coarse grid (~4000 cells) overestimates or underestimates temperatures, and therefore the fine grid is used in the remaining simulations. Among these three resolutions, we also observe that the coarse grid largely underestimates the southwest surface temperature at about 16 h. The reason is that this surface receives less direct solar flux in the coarse grid than in the fine grid as shown on (Fig. 6). In addition, the thermal boundary layer close to the heated surface may need a sufficient resolution to be captured in such a canopy. Hence, we verified that the fine resolution (0.5 m) in the canopy is enough to reproduce a thermal boundary layer thickness of about 2 m.

Fig. 5.

Surface temperature evolution during a diurnal cycle for different mesh resolutions: dashed lines: coarse grid of ~4000 cells; solid lines: fine grid of ~55 000 cells; filled circles: refined grid about 173 000 cells.

Fig. 5.

Surface temperature evolution during a diurnal cycle for different mesh resolutions: dashed lines: coarse grid of ~4000 cells; solid lines: fine grid of ~55 000 cells; filled circles: refined grid about 173 000 cells.

Fig. 6.

Vertical cross section of the direct solar flux in the fluid domain around a building at 16 h for different spatial resolutions: (a) coarse grid and (b) fine grid.

Fig. 6.

Vertical cross section of the direct solar flux in the fluid domain around a building at 16 h for different spatial resolutions: (a) coarse grid and (b) fine grid.

In our model, the DOM was implemented with two angular discretizations: 32 or 128 directions, which influence the prediction of the diffuse solar flux and the infrared flux. In this case, the results obtained with 32 directions are very close to the ones with 128 directions, but the calculation is faster by a factor of 5. So we took 32 directions in the remaining simulations. A time step was introduced for the radiative scheme that is different from the one used for the dynamics. We have tested several radiative time steps: 1, 5, 15, and 30 min and 1 h for a diurnal cycle simulation. The resulting difference between a time step of 1 min and 5 min is small, being less than 1°C. We considered that 5 min was an optimum radiative time step for our simulations. The time step for the dynamics was set to 0.1 s after Milliez and Carissimo (2007). The 24-h simulation in parallel computing on a workstation with eight processors took approximately 4.5 days. We have also tested a dynamical time step of 0.5 s, which reduces the computational time to less than 1 day: the results are close to the ones obtained with a dynamical time step of 0.1 s from 0 to 9 h, but after 10 h an important difference (12%–30%) appears when the wind speed exceeds 6 m s−1 (Fig. 4). In cases in which the wind speed is small, we could set up the dynamical time step to 0.5 s to reduce significantly the CPU time. The full radiative–dynamic coupling remains computationally expensive, in comparison with simpler models, but at this stage the model is intended for research and not for operational applications.

4. Results

a. Simulation of 25 September

Figure 7 shows the evolution of modeled and measured surface temperatures using the force–restore method, with two modeling approaches: 1) radiative only and no convection model (meaning with the convective flux set to zero) and 2) coupled radiative and dynamical model. The diurnal evolutions of the surface temperatures of the top face, southeast face, and northeast face are correctly reproduced by our coupled model. For the northwest face and the southwest face, the simulations show a delay in the morning warming. This delay in warming can be explained by the conduction between the container walls, which contributes to the fast warming in the northwest and southwest faces before they are in the sun but which is neglected in the simulations. In addition, this delay may be inherent to the force–restore method, which overestimates in this case the relaxation to the internal temperature and therefore enhances the thermal inertia. This inertia effect is also observed after sunset for the simulated northwest surface temperature, which shows a delay in cooling. In the afternoon when the atmospheric radiative forcing increases, however, the modeled surface temperatures compare well to the measurements. The comparison with measurements shows a large improvement for the coupled model as compared with the radiative-only model, underlining the importance of accurately including the effect of convection in microscale modeling.

Fig. 7.

Evolution of surface temperature of (a) roof (top), (b) northwest, southeast, (c) northeast, and southwest faces modeled using the force–restore method during a diurnal cycle (crosses: measurements; dashed lines: simulation with radiation only; solid lines: simulation with the dynamic–radiative coupling).

Fig. 7.

Evolution of surface temperature of (a) roof (top), (b) northwest, southeast, (c) northeast, and southwest faces modeled using the force–restore method during a diurnal cycle (crosses: measurements; dashed lines: simulation with radiation only; solid lines: simulation with the dynamic–radiative coupling).

b. Sensitivity to the surface parameters and surface temperature models

The values of the surface parameters were taken from the literature. The values for the thermal properties of metal cover a wide range (Oke 1987), however, and we performed a sensitivity study of the variation of the surface temperatures when varying the parameters in the range given by the literature. Table 2 illustrates that a change in albedo, emissivity, or admittance, in the range given by the literature, can make a difference of about 1°–10°C for the surface temperature. Because the southeast wall is the most exposed to the sun all day, the deviation of the temperature on this wall is the most important.

Table 2.

Deviation of the surface temperatures (°C) during a diurnal cycle by varying albedo, emissivity, and thermal admittance (J m−2 s−0.5 K−1), respectively.

Deviation of the surface temperatures (°C) during a diurnal cycle by varying albedo, emissivity, and thermal admittance (J m−2 s−0.5 K−1), respectively.
Deviation of the surface temperatures (°C) during a diurnal cycle by varying albedo, emissivity, and thermal admittance (J m−2 s−0.5 K−1), respectively.

The results presented above were obtained by using the force–restore method and are now compared with those obtained with the wall thermal model [see section 2b(3)]. The emissivity was chosen after Oke (1987) for corrugated iron. Because some of the containers in the MUST array were painted, we also made a test taking a much higher emissivity (i.e., 0.9). It results in a decrease in surface temperature. This decrease is greater when using the force–restore method, about 10°C for the southeast wall and 5°C for rest of the walls. For the wall thermal model, the influence of higher emissivity on the surface temperature is less significant, the difference being less than 2°C for all of the faces. This is can be explained the fact that the internal temperature has a greater influence than the other surface parameters. The wall thermal model also requires the characteristics of the surface material and the thickness and the thermal conductivity of the wall which were not provided in the data. When choosing a thermal conductivity of 26 W K−1 m−1 for the walls, the resulting conduction remains too high, resulting in an homogenization of the temperature of the five walls (not shown here). To improve the comparisons with the observations, we adjusted the value of the conductivity to 6 W K−1 m−1, which is not that of pure metal, but may be set by assuming an insulating material in the walls. In Fig. 8, we display the evolution of the northwest wall temperature where we used a conductivity of 6 W K−1 m−1 and a thickness of 10 cm for the wall and compare it with the observations and that obtained with the force–restore method. In the morning, the wall thermal model (represented by the circles) is able to simulate accurately the increase in the northwest wall surface temperature at 0600 LT, with no delay, as opposed to the force–restore model (represented by the solid line). An overcooling of the surface temperature appears after sunset, however. This overcooling may be explained by a wrong estimation of the internal temperature by the evolution equation, which is highly dependent on the other computed surface temperatures. Another reason could be an overestimation of the mixing by the turbulence scheme, but we expect this deficiency to be weak. Indeed, the turbulence scheme, which takes into account the stability effects, has been extensively used and was previously validated (Buty et al. 1988; Milliez and Carissimo 2008). The wall thermal model seems more adapted to shipping containers than is the force–restore method, however, and a perspective would be to improve the conduction model by, for instance, implementing a multilayer wall model.

Fig. 8.

The northwest wall temperature evolutions obtained with different surface temperature models (crosses: measurements; filled circles: wall thermal model; solid line: force–restore method).

Fig. 8.

The northwest wall temperature evolutions obtained with different surface temperature models (crosses: measurements; filled circles: wall thermal model; solid line: force–restore method).

5. Discussion: Comparison of three schemes of increasing complexity for predicting surface sensible heat flux

In this section, sensitivity testing is done to compare three schemes used for predicting surface sensible heat flux. The simulated case is based on the MUST geometry with an upstream wind direction of −45°, a reference 10-m wind speed Uref = 4 m s−1, and an initial air temperature of 18°C. The simulation starts at 1200 LT for period of 30 min.

a. Constant hf model

This scheme is usually used in architecture simulation tools (Miguet and Groleau 2002; Asawa et al. 2008). The radiative model in this type of tool is very accurate, usually using a detailed 3D geometry. The convective model is very simplified, however, and the scheme considers a constant transfer coefficient. For the comparison, we take the constant hf as the average value on each wall from the 3D convection model that we presented in section 2b(5). Here, we take hf equal 14.45 W m−2 K−1 for the roof and 6.12 W m−2 K−1 for the walls.

In fact, if instead of taking the same constant hf for all of the walls, we take separate transfer coefficients for each surface (roof 14.45; northwest face 3.94; southeast face 10.78; northeast face 1.38; southwest face 8.35), we can better take into account the orientation of the surface in the wind flow, which decreases the wall surface temperatures by about 2–4 K.

b. One-dimensional hf model (1D hf)

In this model employed in TUF-3D (Krayenhoff and Voogt 2007) and similar to the one used in Masson (2000), the transfer coefficient is calculated based on a simple relationship (Martilli et al. 2002):

 
formula

with u(z) being the vertical wind profile within the canopy. Many authors model this wind profile within the canopy with an exponential law (Cionco 1965; Rotach 1995; Krayenhoff and Voogt 2007). For instance, in TUF-3D, Krayenhoff and Voogt (2007) used an iterative way to find a profile of the exponential form with three coefficients. Here, we model the vertical velocity with the exponential profile of Macdonald (2000), which is well adapted to low-density arrays:

 
formula

where uH is the mean velocity at the top of the obstacles and the constant a is the attenuation coefficient, which is determined by fitting the average wind profile within the obstacle array.

c. Three-dimensional hf model (3D hf)

The full model is three dimensional not only in terms of the radiative exchanges but also the convective exchanges. In this approach, hf is computed by resolving the 3D RANS and energy equations in the whole fluid domain. Coefficient hf is calculated for each subfacet depending on the local friction velocity [Eq. (12)], and the sensible heat flux is calculated with the local air temperature [Eq. (11)].

d. Discussion

Figure 9 illustrates the effect of the three convective schemes by visualizing, successively, the transfer coefficient, the sensible heat flux, and the surface temperature. The three convective schemes show a difference of the sensible flux of approximately 100–240 W m−2 for the southeast face and northeast face. With the constant hf model, the surface temperatures are more homogeneous than in the other two cases. In the MUST configuration, the building array is not dense, and therefore the effects of the shadow and the multireflections are small. That is the reason why the temperatures in the constant hf approach show little difference within each wall. With the 1D hf model, we can obviously see the 1D inhomogeneity of the surface temperatures, which is linked to the exponential law wind profile. The 3D hf model results show the 3D inhomogeneity of the surface temperatures, linked to the inhomogeneity of the 3D wind. On the same face with the same material, we can have a difference of temperature of about 4 K. These results demonstrate the effects of realistically computing the convection fluxes on the surface temperature in the urban areas. Note that in the comparison of the three convective schemes we change only the transfer coefficient and not the air temperature (which is computed for each grid cell of the fluid domain in three dimensions). A simple air temperature model could lead to additional differences.

Fig. 9.

Comparison of three convective models with visualization of the southeast and northeast walls at 1230 LT: (a) transfer coefficient hf (W m−2 K−1), (b) sensible heat flux QH (W m2), (c) surface temperature T (K). Label 1 indicates the constant hf model, label 2 is for the 1D hf model, and label 3 represents the 3D hf model.

Fig. 9.

Comparison of three convective models with visualization of the southeast and northeast walls at 1230 LT: (a) transfer coefficient hf (W m−2 K−1), (b) sensible heat flux QH (W m2), (c) surface temperature T (K). Label 1 indicates the constant hf model, label 2 is for the 1D hf model, and label 3 represents the 3D hf model.

6. Conclusions and perspectives

New atmospheric radiative and thermal schemes were implemented in the atmospheric module of the three-dimensional CFD code known as Code_Saturne. The purpose of this paper was to study in detail the coupling between the radiative scheme and the 3D dynamical model. The model was evaluated with the field measurements from an idealized urban area, the MUST experiment. The coupled model is able to reproduce the evolution of the surface temperatures for different sides of a container within the MUST canopy during a diurnal cycle despite a delay in warming for the northwest and southwest faces at sunrise. The simulations also showed a significant impact of the convective flux on the surface temperatures.

Because the thermal information available in the MUST field is insufficient, sensitivity studies were performed that emphasized the dependence of the model on the parameters describing the building: the properties of the material. In addition, the internal building temperature shows great importance because the buildings are made of metal.

We compared two ways of computing the surface temperature: the force–restore method and a wall thermal model. Because the force–restore method may be more suited for insulated buildings with a near-constant internal temperature (which is not very representative of the MUST containers), we have also tested a one-layer wall model. Using an appropriate evolution equation for the interior buildings, the force–restore shows good results during the afternoon. It induces a delay in warming at sunset, however, and for the northwest wall a delay in cooling, because of the thermal inertia inherent to the method. The one-layer wall thermal model we tested also show some weakness, since it can reproduce the diurnal cycle of the different surface temperature only with a very low thermal conductivity. Nevertheless, it seems more adapted to model shipping containers than is the force–restore method. One proposed improvement is to implement a multilayer 1D thermal model for the walls, as in Masson (2000) or Krayenhoff and Voogt (2007), which may be more adequate for such surfaces. In fact, shipping-container surface temperatures appear difficult to predict with classical urban models, which for instance do not take into account the conduction between the walls, which can be neglected for real buildings but should not be for metal containers, and a 3D conduction wall model may in this case be necessary.

We also compared our 3D modeling approach to estimate the convective exchanges at the surface (which consists of solving the 3D RANS equations in the whole fluid domain) with two modeling approaches found in the literature: the first approach consists in using a constant heat transfer coefficient, and the second approach is to use a 1D equation based on a vertical wind profile within the canopy. The comparisons are made in terms of the convective transfer coefficients, sensible heat fluxes, and surface temperatures. The three schemes give values of the same order of magnitude for the average surface temperature; nevertheless, only the 3D approach can reproduce the inhomogeneous effect on the wind on a surface: the difference of the same wall can reach 4 K.

The simulation of realistic atmospheric conditions in the urban areas made possible by this work can be used for various applications. A first example is to study pollutant dispersal in a low-wind case. A good description of the heat transfer is essential to describe the convective movement of the air in the streets and is very important for air pollution investigations. A second example is the energy balance of buildings. Estimating the convection fluxes in simple models can keep the computing time low and has application in mesoscale studies; nevertheless, at microscale, it can lead to misleading values in the estimation of the energy loss from the buildings to the atmosphere. In this case, a good prediction of the convective flux can be helpful to the management of the energy consumption and a useful tool in building design. The simulations’ results show a significant difference of the parameterizations between taking a simple convection model and describing the physical processes in a realistic way, coupling the 3D dynamics and the radiative processes. The simulations point out the larger difference in surface temperatures at different locations on the same wall.

This study is the first step to validate our dynamic–radiative coupling model. This 3D modeling investigation can bring more detailed information both on radiative and convective fluxes in very local-scale studies. In the MUST case, however, the urban area is idealized. At microscales, small irregularities can break the repeated flow patterns found in a regular array of containers with identical shape. In addition, uncertainties associated with the thickness and the properties of the material of the container wall limit our ability to validate the results. There are also still challenges for modeling in this area. The comparison of three different modeling approaches to estimate the convective exchanges at the surface could be compared with observations if thermal images are available. That is the reason why we will evaluate the coupled dynamic–radiative model on a district of a real urban area with the Canopy and Aerosol Particles Interactions in Toulouse Urban Layer (CAPITOUL) experiment (in the city of Toulouse, France) (Masson et al. 2008; Lagouarde et al. 2010) for which thermal infrared images are available.

Another perspective of this work is to apply the 3D radiative scheme to nontransparent media. Indeed, in many urban applications, the atmosphere between the boundaries can be considered to be transparent and nondiffusive. Nevertheless, when studying smoke dispersion or fog formation and dissipation, absorption and diffusion play an important role. Absorption can already be taken into account by our scheme, and one important perspective of this work is to study radiation in 3D nontransparent media and add the diffusion term in the resolution of the radiative transfer equation.

Acknowledgments

The authors thank the Defense Threat Reduction Agency (DTRA) for providing access to the MUST data.

REFERENCES

REFERENCES
Archambeau
,
F.
,
N.
Méchitoua
, and
M.
Sakiz
,
2003
:
Code_Saturne: A finite volume code for the computation of turbulent incompressible flows—Industrial applications
.
Int. J. Finite Vol.
,
1
,
1
62
.
Arnfield
,
A. J.
,
2003
:
Two decades of urban climate research: A review of turbulence, exchanges of energy and water, and the urban heat island
.
Int. J. Climatol.
,
23
,
1
26
.
Asawa
,
T.
,
A.
Hoyano
, and
K.
Nakaohkubo
,
2008
:
Thermal design tool for outdoor space based on a numerical simulation system using 3D-CAD
.
Build. Environ.
,
43
,
2112
2123
.
Biltoft
,
C. A.
,
2001
:
Customer report for Mock Urban Setting Test
.
West Desert Test Center, U.S. Army Dugway Proving Ground, DPG Doc. WDTC-FR-01-121, 58 pp
.
Brook
,
D. R.
,
N. V.
Beck
,
C. M.
Clem
,
D. C.
Strickland
,
I. H.
Griffits
,
D. J.
Hall
,
R. D.
Kingdon
, and
J. M.
Hargrave
,
2002
:
Validation of the Urban Dispersion Model (UDM)
.
Proc. Eighth Int. Conf. on Harmonisation within Atmospheric Dispersion Modelling for Regulatory Purposes, Sofia, Bulgaria, National Institute of Meteorology and Hydrology of Bulgaria, 8–12
.
Buty
,
D.
,
J. Y.
Caneill
, and
B.
Carissimo
,
1988
:
Simulation numerique de la couche limite atmospherique en terrain complexe au moyen d’un modele mesometeorologique non hydrostatique: Le code mercure (Numerical simulation of the atmospheric boundary layer in complex terrain by means of a nonhydrostatic model: The Mercury code)
.
J. Theor. Appl. Mech.
,
7
,
35
62
.
Camelli
,
F. E.
,
R.
Lohner
, and
S. R.
Hanna
,
2005
:
VLES study of MUST experiment
.
Proc. 43rd AIAA Aerospace Sciences Meeting and Exhibit, Reno, Nevada, American Institute of Aeronautics and Astronautics, 1–19
.
Cionco
,
R. M.
,
1965
:
A mathematical model for air flow in a vegetative canopy
.
J. Appl. Meteor.
,
4
,
517
522
.
Deardorf
,
J. W.
,
1978
:
Efficient prediction of ground surface temperature and moisture with inclusion of a layer of vegetation
.
J. Geophys. Res.
,
83
,
1889
1903
.
Dupont
,
S.
, and
P.
Mestayer
,
2006
:
Parameterization of the urban energy budget with the submesoscale soil model
.
J. Appl. Meteor. Climatol.
,
45
,
1744
1765
.
Eichhorn
,
J.
, and
M.
Balczo
,
2008
:
Flow and dispersal simulations of the Mock Urban Setting Test
.
Proc. 12th Int. Conf. on Harmonisation within Atmospheric Dispersion Modelling for Regulatory Purposes, Cavtat, Croatia, Meteorological and Hydrological Service of Croatia, 67–72
.
Fiveland
,
W. A.
,
1984
:
Discrete-ordinates solutions of the radiative transport equation for rectangular enclosure
.
J. Heat Transfer
,
106
,
699
706
.
Garratt
,
J. R.
,
1992
:
The Atmospheric Boundary Layer
.
Cambridge University Press, 316 pp
.
Gastellu-Etchegorry
,
J.
,
2008
:
3D modeling of satellite spectral images, radiation budget and energy budget of urban landscapes
.
Meteor. Atmos. Phys.
,
102
,
187
207
.
Gastellu-Etchegorry
,
J.
,
E.
Martin
, and
F.
Gascon
,
2004
:
DART: A 3-D model for simulating satellite images and surface radiation budget
.
Int. J. Remote Sens.
,
25
,
75
96
.
Grimmond
,
C. S. B.
, and
T. R.
Oke
,
1999
:
Heat storage in urban areas: Local-scale observations and evaluation of a simple model
.
J. Appl. Meteor.
,
38
,
922
940
.
Hanna
,
S. R.
,
S.
Tehranian
,
B.
Carissimo
,
R. W.
Macdonald
, and
R.
Lohner
,
2002
:
Comparisons of model simulations with observations of mean flow and turbulence within simple obstacle arrays
.
Atmos. Environ.
,
36
,
5067
5079
.
Johnson
,
G. T.
,
T. R.
Oke
,
T. J.
Lyons
,
D. G.
Steyn
,
I. D.
Watson
, and
J. A.
Voogt
,
1991
:
Simulation of surface urban heat islands under ‘ideal’ conditions at night. Part I: Theory and tests against field data
.
Bound.-Layer Meteor.
,
56
,
275
294
.
Kanda
,
M.
,
T.
Kawai
,
M.
Kanega
,
R.
Moriwaki
,
K.
Narita
, and
A.
Hagishima
,
2005
:
A simple energy balance model for regular building arrays
.
Bound.-Layer Meteor.
,
116
,
423
443
.
Krayenhoff
,
E. S.
, and
J. A.
Voogt
,
2007
:
A microscale three-dimensional urban energy balance model for studying surface temperatures
.
Bound.-Layer Meteor.
,
123
,
433
461
.
Lagouarde
,
J. P.
,
A.
Hénon
,
B.
Kurz
,
P.
Moreau
,
M.
Irvine
,
J.
Voogt
, and
P.
Mestayer
,
2010
:
Modelling daytime thermal infrared directional anisotropy over Toulouse city centre
.
Remote Sens. Environ.
,
114
,
87
105
.
Lemonsu
,
A.
,
C. S. B.
Grimmond
, and
V.
Masson
,
2004
:
Modelling the surface energy balance of an old Mediterranean city core
.
J. Appl. Meteor.
,
43
,
312
327
.
Liu
,
J.
,
H. M.
Shang
, and
Y. S.
Chen
,
2000
:
Development of an unstructured radiation model applicable for two-dimensional planar, axisymmetric, and three-dimensional geometries
.
J. Heat Transfer
,
66
,
17
33
.
Louis
,
J.
,
1979
:
A parametric model of vertical eddy fluxes in the atmosphere
.
Bound.-Layer Meteor.
,
17
,
187
202
.
Macdonald
,
R. W.
,
2000
:
Modelling the mean velocity profile in the urban canopy layer
.
Bound.-Layer Meteor.
,
197
,
25
45
.
Martilli
,
A.
,
A.
Clappier
, and
M. W.
Rotach
,
2002
:
An urban surface exchange parameterisation for mesoscale models
.
Bound.-Layer Meteor.
,
104
,
261
304
.
Masson
,
V.
,
2000
:
A physically based scheme for the urban energy budget in atmospheric models
.
Bound.-Layer Meteor.
,
94
,
357
397
.
Masson
,
V.
,
C. S. B.
Grimmond
, and
T. R.
Oke
,
2002
:
Evaluation of the Town Energy Balance (TEB) scheme with direct measurements from dry districts in two cities
.
J. Appl. Meteor.
,
41
,
1011
1026
.
Masson
,
V.
, and
Coauthors
,
2008
:
The Canopy and Aerosol Particles Interactions in Toulouse Urban Layer (CAPITOUL) experiment
.
Meteor. Atmos. Phys.
,
102
,
135
157
.
Miguet
,
F.
, and
D.
Groleau
,
2002
:
A daylight simulation tool for urban and architectural spaces—Application to transmitted direct and diffuse light through glazing
.
Build. Environ.
,
37
,
833
843
.
Milliez
,
M.
,
2006
:
Modélisation micro-météorologique en milieu urbain: Dispersion des polluants et prise en compte des effets radiatifs (Micrometeorological modeling in urban environment. Dispersion of pollutants and the taking into account of radiative effects). Ph.D. thesis, École des Ponts ParisTech, 228 pp. [Available online at http://cerea.enpc.fr/fich/these_milliez.pdf.]
Milliez
,
M.
, and
B.
Carissimo
,
2007
:
Numerical simulations of pollutant dispersion in an idealized urban area, for different meteorological conditions
.
Bound.-Layer Meteor.
,
122
,
321
342
.
Milliez
,
M.
, and
B.
Carissimo
,
2008
:
CFD modelling of concentration fluctuations in an idealized urban area, for different meteorological conditions
.
Bound.-Layer Meteor.
,
127
,
241
259
.
Milliez
,
M.
,
L.
Musson-Genon
, and
B.
Carissimo
,
2006
:
Validation of a radiative scheme for CFD modelling of heat transfers between buildings and flow in urban canopies
.
Preprints, Sixth Int. Conf. on Urban Climate, Göteborg, Sweden, International Association of Urban Climate, 266–269
.
Mills
,
G.
,
1996
:
An urban canopy-layer climate model
.
Theor. Appl. Climatol.
,
1997
,
229
244
.
Oke
,
T. R.
,
1987
:
Boundary Layer Climates
.
2nd ed. Routledge, 435 pp
.
Rotach
,
M. W.
,
1995
:
Profiles of turbulence statistics in and above an urban street canyon
.
Atmos. Environ.
,
29
,
1473
1486
.
Truelove
,
J. S.
,
1987
:
Discrete-ordinates solutions of the radiative transport equation
.
J. Heat Transfer
,
109
,
1048
1051
.
Yee
,
E.
, and
C. A.
Biltoft
,
2004
:
Concentration fluctuations measurements in a plume dispersing through a regular array of obstacles
.
Bound.-Layer Meteor.
,
111
,
363
415
.