Use of wind pressure coefficients to simulate natural ventilation and building energy for isolated and surrounded buildings

.


Introduction
Wind pressure coefficients (C p ) are key inputs for natural ventilation calculation in building energy simulations (BES) and multi-zone airflow models (e.g.AirflowNetwork (AFN) in EnergyPlus (EnergyPlus-AFN) [1], MacroFlo in IES-VE [2], CONTAM and COMIS linkages with TRNSYS [3]).C p is the nondimensional ratio of wind pressure on the building surface to the dynamic pressure in the upstream undisturbed flow [4] but is defined differently depending on the height of free stream dynamic pressure [5]: where H is a reference height (m, building roof/eave height), or alternatively the local opening height z (m) is used [5]: where P w (z) is the wind pressure (Pa) measured on the building facet at height z (m), U free is the wind speed (m s − 1 ) in the upstream undisturbed flow at H or z, and ρ is the outdoor air density (kg m − 3 ) which is assumed constant.Typically, C p is calculated as the average value across the entire building facet facing the flow (i.e.surface mean).C p is widely applied in studies of natural ventilation potential [6][7][8][9], cooling energy savings [10][11][12][13], indoor thermal comfort and overheating [14][15][16][17], and other applications like the solar chimney [18] and windcatcher [19].Commonly the C p data sources used in BES (Table 1) are from Ref. [21]: primary sources (e.g., full-scale experiments, wind tunnel experiments and CFD simulations for a specific building of interest); and secondary sources (e.g.databases with generic building archetypes derived often from wind tunnel experiments).In databases (Table 1), C pr rather than C pl data are provided and are the default values used in BES.
C p data, together with the reference wind conditions, are regarded as the major sources of uncertainties in multi-zone airflow models [26,27].Therefore, applying C p data correctly is important for more accurate estimation of natural ventilation rates.C p data are dependent on many factors including the height of free stream dynamic pressure measured and the vertical wind profile which is modified by the surrounding buildings [22,24].However, these dependencies are often overlooked in BES, which could cause biases.Here we focus on EnergyPlus, as it is one of the most widely used open-source BES tools.Our aim is to critically explore potential biases of EnergyPlus-AFN simulation in three comparative C p application scenarios (Fig. 1).
Scenario 1: In EnergyPlus the surface averaged C pr (C pr ) data are usually used with wind speed at the opening height U free (z).This is inconsistent (Eq.( 1) cf.Eq. ( 2)) and can cause biases since U free (z) should be used with the surface averaged C pl (C pl ) (Fig. 1a).EnergyPlus options allow the use of either provided [5,24] or user-supplied C p data.In the latter case, the user also needs to indicate the height the C p values are for (i.e., opening or reference), with the EnergyPlus default (i.e., if not modified) being the opening height [28].If the provided default C p values are used, the opening height will be used (not explicitly stated in Ref. [28] but found in the source code [29] and mentioned by Refs.[13,30,31]).Hence, the provided C pr value is used with wind speed at the opening height (z) instead of the reference roof height (H).This inconsistency between C pr and freestream wind speed will cause biases in the wind pressure calculations if not corrected.Given this, we review current studies for choosing C p if modelling natural ventilation in buildings using EnergyPlus-AFN (Table 2).Typically, if the detailed C p settings are not mentioned, the supplied default C pr values are assumed to be used.

A
Opening area (m  Of the studies stating the C p settings (fully or partly), only four [31,33,40,43] of 27 (Table 2) use the correct combination of C p and free stream wind speed (i.e.C pr with U free (H) or C pl with U free (z)).The other studies (where both C p and U free are clear) all use C pr with U free (z).Thus, this common bias (23 of 27 studies) caused by using C pr with U free (z) needs to be assessed.Scenario 2: C pr is defined using pressure and wind profiles from wind tunnel or CFD studies.If the wind profile in building energy models have systematic differences compared to the wind profile used to derive the pressure coefficients, then using the unmodified C pr data will cause systematic errors in predicted pressure values (Fig. 1b).C pl is calculated with the wind speed at the same height as the wind pressure, so is not a function of the wind profile [4,24], and is directly applicable in BES with a different wind profile from the one used to derive C pl in the wind tunnel or CFD study.However, this does not apply to C pr because it is based on the wind speed at a reference height [55,56].For example, it is possible for different wind profiles to have the same wind speed at the reference height H. Neglecting differences in vertical wind profiles can cause biases.Potentially, this is a large problem as C pr is widely used (Tables 1 and 2), and wind profiles in BES are normally different from the C p source experiment wind profiles especially when using the secondary sources (11 of 27 studies in Table 2).When the wind profiles for both the C pr source and the BES are known, the C pr data can be modified appropriately.
Scenario 3: If BES are combined with urbanised wind speed from urban canopy models, the free stream wind C p values should be also accounted for the influence of the surrounding buildings (Fig. 1c).With increasing attention on urbanization and the impact of urban climate on building performance, efforts have been made to integrate BES with urban land surface or canopy models, such as combining EnergyPlus with Surface Urban Energy and Water Balance Scheme (SUEWS) [57][58][59] and the Vertical City Weather Generator (VCWG) [60].Urban canopy models can modify meteorological variables to better account for the impact of buildings on local climate which will influence wind pressure calculation.For example, SUEWS provides a mean neighbourhood vertical profile of wind speed [58,61].This differs from the undisturbed wind used in Eqs. ( 1) and ( 2) which are default weather data inputs in EnergyPlus.In this case, C p values need to be corrected if the local wind speed is used (Fig. 1c).This has been largely overlooked in existing urban-building coupling energy simulation studies (e.g., use of C p for surrounded case and disturbed local wind speed as reviewed by Ref. [62]).
The objectives of this study are to: (1) quantify the bias arising from using inconsistent reference height and wind pressure coefficient combinations; (2) assess the bias arising from the inconsistency of approaching wind profiles between BES and the source deriving the wind pressure coefficients; and (3) discuss the correction of wind pressure coefficients when combining building energy simulation tools and urban climate models.

Methods
To analyse the use of wind pressure coefficients in building energy simulation (BES), we use the BES tool EnergyPlus v9. 4 [63] with Airflow Network (AFN) for ventilation rate calculation, and urban land surface model SUEWS [57,64] to simulate urban wind profiles.The wind pressure coefficient (C p ) data are obtained from the Tokyo Polytechnic University's aerodynamic database for low-rise buildings [25] under isolated and surrounded scenarios.
Each floor of the building has a 3 m × 2 m window on the north and south facing walls, of which the upper 1/3 area hinged and openable to 20 • for cross ventilation.For consistency with weather data, the building envelope thermal characteristics are set using the current local Shanghai building code [67].This has overall heat transfer coefficients (U-values) of 0.39 W m − 2 K − 1 for the roof, 0.54 W m − 2 K − 1 for the external wall, 0.46 W m − 2 K − 1 for the floor, and 1.77 W m − 2 K − 1 for the windows.All windows are assumed to have 15% openable area and a discharge coefficient (C d ) of 0.61.
We calculate the natural ventilation rate, indoor overheating risk and energy saving potential.For the naturally ventilated mode, all windows are always open.The overheating risk is assessed using the Category II Chinese adaptive thermal model comfort corresponding to 75% satisfaction [68].For the southern zone (i.e., applicable for Shanghai) the upper (T max ) and lower temperature limits (T min ) are [68]: where the running mean outdoor temperature T rm is: where k is a constant between 0 and 1, with 0.8 used as recommended [69], and T od-n is the daily mean outdoor temperature n days ago ( • C).
The air-conditioning mode follows the Code for Thermal Design of Civil Building recommendation of heating (18 • C) and cooling (26 • C) setpoints [70].Windows can be open if air-conditioning is off and outdoor temperature is lower than indoor temperature.The cooling energy saving is calculated as the difference between the energy demand in hybrid mode (natural ventilation together with air conditioning) and fully air-conditioned mode.
The indoor overheating metrics are hours and degree hours exceeding the upper limits of temperatures (T max ) [71].Local outdoor weather data required for EnergyPlus simulations for the two neighbourhoods are generated using SUEWS with the vertical profiles option of air temperature, relative humidity, and wind speed [58].Air temperature and wind profiles evaluations using observations at three sites have reasonable accuracy [58,61].For the neighbourhood cases, solar shading and inter-building external longwave radiative exchanges are also considered [72].
The normalised mean bias error (nMBE) assessment metric is used to compare two cases (x; y) for each of the three scenarios (Section 1, Fig. 1): where x i is the results from a consistent combination of C p and level of U, and y i is the results from the inconsistent combination of these for each hour i in the year (total of N = 8760 h).The ASHRAE-14 Guideline [73] acceptable uncertainty limits for building energy simulation programmes is the nMBE needs to be within ±10% for hourly data.

C pl & U free (z) vs C pr & U free (z) for an isolated building
First, we need to derive a relation between surface averaged C pl (C pl ) and C pr (C pr ).To obtain the surface averaged C p over the facet one can either calculate C p at several locations and average, or define C pl as the constant value that gives the correct total wind pressure over the facet [4].The latter is analogous to Santiago and Martilli's [74] approach used for vertically distributed drag modelling of urban canopies, where the surface-averaged drag coefficient is the value giving the correct total drag over the building facet.Akins and Cermak [4] suggested the differences between these two techniques are minimal, so we assume surface averaged C p calculated with both methods are the same.We take the correct total wind pressure approach, which is defined by rearranging Eq. (1) or Eq 2 for P w (z) and integrating over the facet.For the pressure coefficient defined with velocity at the local opening height (Eq.( 3)) one finds: and for the pressure coefficient defined with velocity at the building height (Eq.( 1)) one finds: By only integrating over height it has been assumed that C p variations in the horizontal can be neglected or that P w (z) has first been horizontally averaged across the facet.The approach is practical since velocity profiles in building energy models normally have vertical variation, so horizontal variation of C p is not included.
To get the same average wind pressure over the building facet 1 H ∫ H 0 P w (z)dz, one can combine Eqs. ( 6) and ( 7) so that: which after rearranging becomes: To find C pl , values for C pr and an equation for the wind speed are required.C pr data from the TPU [25] database are used.Commonly in wind tunnel experiments (e.g.Ref. [25], a power law is used to describe the undisturbed wind speed at height z): where U ref is the reference wind speed at height z ref defined within the experiment, and the exponent α is an empirically derived coefficient.
The vertically averaged wind speed is given by: Substituting Eqs. ( 10) and (11) into Eq.( 9): Using Eq. ( 12) with C pr data from the TPU database (vertical profiles of original C pr data are shown in Supplementary Material Fig. SM.1a), the results for an isolated reference building are given in Table 3, as are the C pl data.EnergyPlus is used with the power law velocity profile that is the same as in the wind tunnel experiment in the TPU database assuming a suburban terrain and an exponent α of 0.2.EnergyPlus simulations are conducted for Shanghai in 2018 under naturally ventilated and mechanical cooling/heating modes as described in Section 2. The normalised mean bias errors (nMBE) for the ventilation rate, indoor overheating risks and cooling/heating energy demand are calculated between C pl & U free (z) and C pr & U free (z).
The C pr with U free (z) underpredicts the annual air change per hour (ACH) ventilation rate (Fig. 3), with a nMBE of − 15.5% (i.e., exceeding the ASHRAE-14 acceptable limit of ±10%).The annual overheating hours and degree hours are overpredicted by 11.9% and 12.9%, respectively.Such differences in overheating are found to largest during three consecutive heatwave days (13-15 July 2018, Fig. 4).During this period the C pr & U free (z) case overpredicts the overheating hours and degree hours by 19.5% and 12.4%, respectively.Given the smaller ventilation rate, the nMBE in cooling energy saving is − 10.5%, again exceeding the ASHRAE-14 limit.This suggests confusing C pl with C pr should be avoided when modelling ventilation rates and indoor overheating risks of naturally ventilated buildings, as well as the resultant cooling energy saving.

C pr with wind profiles: wind tunnel vs outdoor
Swami and Chandra [24] and Akins et al. [5] suggest C pl is not a function of the wind profile given it is based on the wind at the opening height z, which is the same as pressure measurement height.After further testing (section SM.2), the results suggest that C pl independence on wind profile exponent α is acceptable.
However, C pr obviously depends on the wind profile.Therefore, when the wind profile in the building energy simulation is different from the wind tunnel experiment where the C p data are derived, C pl can be used directly without further corrections.If only undisturbed wind speed at height H is available (e.g., TMY (typical meteorological year) wind speed data at 10 m), the C pr values should be corrected.The EnergyPlus outdoor wind profile module determines the approaching wind speed profile U(z) as [20]: It is calculated with the wind speed measured at a meteorological station U met (i.e., weather data input).Standard World Meteorological Organisation (WMO) wind speed measurement height (z met ) is 10 m above ground level [75].δ refers to the height where the vertical gradient of wind speed is assumed to become constant [58].Typical values are given in Ref. [20] for different terrain types.For an isolated building in open country terrain α = α met = 0.14 and δ = δ met = 270 m.Following Eq. ( 8), assuming C pl data from the TPU database can be used directly (i.e.C pl,WT = C pl,EP ) the average wind pressure on the building facet (Fig. 4) is calculated as and upon rearranging one finds And similar to Eq. ( 12): where C pr,EP and C pl,EP are used with the EnergyPlus wind speed (U EP ), and C pr,WT and C pl,WT are obtained from the wind tunnel experiment.
To quantify the bias from neglecting the impact of vertical wind profile on C pr , we model the isolated reference building in EnergyPlus, and calculate the ventilation rate, indoor overheating risks and cooling energy demand with both modified C pr,EP and unmodified C pr,WT obtained from the TPU database directly (Table 4).In the TPU database the wind profile exponent is α WT = 0.2.The default EnergyPlus wind profiles exponents (α EP ) are 0.14, 0.22 and 0.33, for open, rough and urban terrain [20], respectively.
Results (Fig. 5) show there are biases if an unmodified C pr,WT is used along with varying α EP .When α EP < α WT , C pr,WT is smaller than C pr,EP , hence the ventilation rate is underpredicted when using C pr,WT .At α EP = 0.14, the nMBE in ACH is − 4.5%, which is within the acceptable range of ASHRAE-14.The annual indoor overheating hours and degree hours are overpredicted by 2.3% and 3.7%, respectively.During the heatwave (13-15 July 2018) the overpredictions slightly increase to 8.0% and 4.4%.nMBE in cooling energy saving is − 3.3%.With α EP = 0.22 and 0.33, the ventilation rates are overpredicted resulting in underpredicted overheating risks and overpredicted cooling energy saving (Fig. 5).Although all the biases are smaller than the ±10% ASHRAE-14 limit, they would increase in rougher terrain as the wind profile exponent increases.

Wind pressure coefficients for urban climate models
When local outdoor weather data are derived from urban weather/ climate models, the influence of the neighbourhood buildings are considered, but when calculating building facet wind pressure an 'undisturbed' flow is assumed in Eqs. ( 1) and (2).For example in SUEWS, the horizontally averaged neighbourhood wind speed is calculated for the roughness sublayer (RSL) based on a modified MOST (Monin-Obukhov similarity theory) approach [61].Similar methods are used in other models like the Vertical City Weather Generator [60].The

Table 3
Surface-averaged wind pressure coefficients based on the reference height H (C pr ) and opening height z (C pl ) by wind angle related to the facet (0 • is when wind is normal to the facet).C pr is obtained from the TPU [25] database.C pl is calculated by substituting C pr into Eq.( 12).advantages of using RSL wind are in calculating convective heat transfer on building surfaces and single-sided ventilation.But the cross-ventilation calculation will be biased if the urbanised wind speed is not used with corrected C p values.
In SUEWS, the undisturbed wind profile can be modelled with an open rural setting (U r ), and the disturbed RSL wind profile can be modelled using different urban settings (U u ).To get the average wind pressure across the building facade: where C pl,u is the corrected wind pressure coefficient for use alongside with the RSL wind (U u ) assuming the vertical profiles of U u and U r have the power law format of Eq. ( 13).Details of obtaining power law vertical profiles of U u and U r applicable for use with EnergyPlus are given in Ref. [58].Eq. ( 19) can be re-written as: Since the RSL wind have lower velocities than the undisturbed wind, C pl,WT values need to be scaled to larger magnitudes (C pl,u ) to obtain the same wind pressure.To quantify the biases of using the disturbed RSL wind speed U u with C pl,WT , we consider two idealised neighbourhoods with aligned buildings in EnergyPlus with plan area fractions of λ P = 0.3 and 0.6, with the surface-averaged C p values given in Table 5  For the neighbourhood with a λ P of 0.3, using C pl,WT with the RSL wind U u (z) largely underpredicts the ventilation rate, with an annual nMBE of − 19.0%, exceeding the ASHRAE-14 acceptable uncertainty limits.The annual indoor overheating hours and degree hours are overpredicted by 5.9% and 13.2%, respectively.During the heatwave

Table 4
Surface-averaged wind pressure coefficients from the wind tunnel experiment (C pr,WT ) [25] and corrected with Eq. ( 17) (C pr,EP ).The wind angles are relative to the surface (0 • refers to wind blowing perpendicular to the facet).(13-15 July 2018), the indoor overheating hours difference is 0, because in both cases the indoor air temperature exceeds the maximum temperature threshold throughout (Fig. 6a), but the overprediction in degree hours increases to 18.4%.The nMBE in cooling energy saving is − 14.0%.Additionally, the ventilation rate is simulated with C pl,WT and the free stream wind U r (z) to evaluate Eq. (18).Results (Fig. 7), as expected, suggest that using C pl,WT with U r (z) and C pl,u with U u (z) give very similar results (nMBE = 0.8%).The small differences are possibly due to the buoyancy-driven ventilation dominating, since the wind speed input to EnergyPlus is also used for calculations of both convection and indoor air temperature.
Generally for the λ P = 0.6 neighbourhood, the biases are slightly smaller (cf.λ P = 0.3) because of the lower wind speeds.The nMBE for the ventilation rate is − 16.2%, while the annual indoor overheating hours and degree hours are overpredicted by 4.7% and 9.0%, respectively.However, during the heatwave (Fig. 6b) when ventilation rates are small, differences in overheating hours become 20.3%.The nMBE in cooling energy saving is − 9.6%.
In summary, when modelling naturally ventilated buildings using urbanised wind speeds, correcting the C p data correspondingly is important.Increasing λ P can lead to slightly smaller annual biases, but during heatwaves when the natural ventilation rates are small, larger biases can be seen.

Discussion
There are various assumptions and approximations made for C p to simplify the calculation of wind pressure on building facets that can cause various uncertainties.Some of these have been assessed previously, such as those linked to the surface averaged values and different data sources [76,77].Most common sources of C p data provide surface averaged values based on the reference height (C pr ), rather than being based on the opening height (C pl ).Given the definition (Eq.( 7)), C pr data needs to be corrected in some circumstances, but this appears to have been overlooked in most existing studies.In this study we explore three scenarios to quantify biases from inconsistent combination of C p value and wind speed.
In each scenario we find critical differences, which impact the resulting predictions especially of ventilation rates and indoor overheating risks for naturally ventilated buildings.These findings confirm that natural ventilation rate calculations are sensitive to the C p and wind data used.Notably, we revise the relation between the C pr and C pl that has often been neglected in building energy simulations.Our results demonstrate the importance of modifying C p data for wind conditions, including the wind speed height, wind profile and terrain surface type (e.g.extensive grass -'undisturbed', neighbourhoods at different λ P -'disturbed' or 'urbanised').
There are limitations in our work.Although surface averaged C p data are widely used, their errors (cf.local C p data) are assumed to be relatively smaller if openings are located in the facet centre instead of edges where extreme values occur [76].Hence, we only consider windows located in the centre of each facet.We consider only one climate type, but expect that relative results should be similar across different

Fig. 1 .
Fig. 1.Three scenarios (Section 1) to determine wind pressure coefficients (C p ) all assume they are surfaceaveraged (C pl or C pr ) but with different wind profiles (U): (a) Scenario 1 are calculated with free stream wind speed at the opening height (U free (z)), (b) Scenario 2: C pr is derived from wind tunnel (WT) wind profile (C pr,WT ) or from the EnergyPlus (EP) building energy simulation (C pr,EP ); and (c) Scenario 3: C pl is based on disturbed urbanised (u) wind speed (C pl,u ) rather than an undisturbed free stream wind speed (C pl,WT ).

Fig. 3 .
Fig. 3. Distribution of annual ACH (air change per hour, N = 8760) calculated with surface averaged wind pressure coefficients based on the opening height z (C pl ) and reference height H (C pr ), with interquartile range (box), median (orange line) and 5th and 95th percentiles (whiskers).(For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.) (original C pr,WT vertical profiles are shown in Fig. SM.1b, c).

Fig. 4 .
Fig. 4. Diurnal changes during three-day heatwave (13-15 July 2018) period in Shanghai for the upper floor (Fig. 2a) (a) indoor air temperatures, and (b) ventilation rates in ACH when calculated using the C pl and C pr .i.e., C pl & U free (z) vs C pr & U free (z).

Fig. 5 .
Fig. 5. Normalised mean bias errors (nMBE) linked to using unmodified C pr,WT from the wind tunnel experiment compared to modified C pr,EP for three EnergyPlus wind profiles exponents (α EP ) for (a) ventilation rate in ACH, (b) difference in annual overheating risks and (c) cooling energy saving.

Fig. 6 .
Fig. 6.As Fig. 4, but using surface-averaged wind pressure coefficients from wind tunnel experiments (C pl,WT ) [25] and corrected with Eq. (20) (C pl,u ) at for two neighbourhoods when λ P is (a, c) 0.3 and (b, d) 0.6, for (a, b) air temperature and (c, d) ventilation rate in ACH.

Fig. 7 .
Fig. 7.As Fig. 3, but for modelling ventilation rate (as ACH) with three combinations of C pl and U(z).C p,WT & U r (z) and C p,u & U u (z) give the similar results but vary slightly due to differences in buoyancy-driven ventilation, whereas assuming C p,WT & U u (z) is inconsistent and therefore biased.

Table 1
Summary of commonly used databases (DB#) of wind pressure coefficients (C [21]in building energy simulations (BES) either using reference height (C pr ) or the local height (C pl ).C pr and C pl refer to surface averaged C pr and C pl data, respectively.*Cpldataarealsocalculatedbut not used in ASHRAE Handbook of Fundamentals[20]and EnergyPlus[1].In places information is not given (NG).Examples of where the data source are used in a BES tool as a default are given.Modified after Cóstola et al.[21].

Table 2
[21]ary of types of C p data and wind speed height used in EnergyPlus and Airflow Network (AFN) studies.All C p values are surface averaged values.Following[21], sources are either primary (1 • ) from wind tunnel (WT) experiments and computational fluid dynamics (CFD) simulations, or secondary (2 • ) from published databases (DB, with references to Table1indicated by #) or analytical tools for generic building archetypes.EnergyPlus-AFN studies that do not indicate the C p method used, are assumed to use pre-provisioned C pr values with wind speed at the opening height z.H refers to the reference building height.Sometimes information is not given (NG).
X.Xie et al.