不同建筑环境下平均风速剖面的数值模拟

摘要

基于计算流体力学（CFD）方法对具有不同建筑环境参数的城市上空的风场进行系统性的数值模拟，并对平均风速剖面指数的变化规律进行讨论。研究结果表明，建筑平均高度对平均风速剖面指数的影响较大，建筑密度对平均风速剖面指数的影响较小。文中还给出了以建筑平均高度和建筑密度为自变量的平均风速剖面指数经验公式。

关键词

平均风速剖面指数; 建筑环境参数; 城市风场; 数值模拟

近40年，高层建筑、高耸结构等风敏感结构在我国城市大量涌现，使得地表对大气流动的摩擦阻力大大增加，大气边界层风剖面的形状发生显著改变。为了研究边界层的形成机理，一些平均风速剖面的研究通过现场实测^［1-3］和风洞试验^［4-7］的方式开展。除了实测和试验之外，城市风场的数值模拟方法也逐渐发展成熟^［8-12］。城市风场的研究涉及城市室内外环境中热量和物质的传递，对于建筑舒适度以及污染物扩散等研究有着重要的意义。因此，有关城市风场的数值模拟关注点一般在于城市街区内部某处的风速或者建筑表面风压，而对于城市高空风场的演变规律缺少相应的研究。

相对国外城市而言，我国城市普遍呈现出建筑物高度更高更密集的特征，但我国的建筑结构荷载规范（GB50009—2012）^［13］对城市上空平均风速剖面指数的确定主要参考发达国家的相关规定，与我国密集城市的实际风速剖面指数存在较大的差异。另外，在实际工程中，工程人员需要对地貌类别进行准确地判断，进而对风敏感结构的风荷载进行计算，但各国相关规范对地貌类别的描述都是定性的，没有给出定量的确定方法，导致工程人员对地貌类别的确定存在很大的主观性，也可能导致工程设计过分保守或者不够安全的情况出现。本文对具有不同建筑环境参数的城市上空的风剖面进行系统性的数值模拟，给出通过建筑环境参数判断城市上空平均风速剖面指数的公式，为解决实际工程问题提供参考。

1 模拟方法验证

为了验证城市风场数值模拟的准确性，本文采用Brown等^［6］在美国国家环境保护局气象学风洞中得出的试验数据进行对比。该风洞试验段尺寸18.3 m×3.7 m×2.1 m，试验采用脉冲热线风速计采集风速时程。粗糙元阵列采用7×11的立方体，立方体规格为0.15 m×0.15 m×0.15 m，粗糙块之间的间隔也为0.15 m。风洞试验的入口采用指数律风速剖面，平均风速剖面指数为0.16。

图1为计算域以及边界条件参数设置，根据风洞试验情况，计算域大小设置为4.95 m×1.50 m×1.20 m。在图中A和B两个位置处设置观测剖面，与风洞试验相同位置处的试验结果进行对比。

图1  计算域及边界条件

Fig.1  Computational domain and boundary condition

尽管四面体网格在收敛性和计算精度上要弱于六面体结构化网格，为了适应城市建筑造型，城市风场模拟中仍常采用四面体网格进行划分。所以本文在建筑周边较小区域内布置四面体网格来适应建筑外形，而在外部区域使用六面体网格减少网格数量。最小网格长度设置为0.01 m，网格总量达到了286万（图2）。

图2  网格划分情况

Fig.2  Scheme of grid discrete

表1为该算例的边界条件设置。该算例采用ANSYS Fluent18.0模拟平台进行计算，湍流模型采用Standard k-ε模型，求解设置压力速度耦合方式为SIMPLEC，动量方程和湍流模型方程非线性对流项采用second order upwind离散格式，压力插值格式采用second order，梯度插值方法采用least squares cell based。所有变量和连续性方程的残差收敛标准设置为10^-6。数值模拟入口采用风洞试验数据拟合：$u_{\mathrm{i}\mathrm{n}\mathrm{l}\mathrm{e}\mathrm{t}}=\mathrm{3}\mathrm{m}·{\mathrm{s}}^{-\mathrm{1}}{(z/\mathrm{0.15}\mathrm{m})}^{\mathrm{0.16}}$；$k=\mathrm{0.2}{\mathrm{m}}^{\mathrm{2}}·{\mathrm{s}}^{-\mathrm{2}}$；$\epsilon ={C_{\mu }}^{\mathrm{3}/\mathrm{4}}{k}^{\mathrm{3}/\mathrm{2}}/\kappa z$。其中，$u_{\mathrm{i}\mathrm{n}\mathrm{l}\mathrm{e}\mathrm{t}}$、$k$、$\epsilon$分别为风速、湍动能和湍流耗散率；$z$为高度；$C_{\mu }$、$\kappa$为湍流模型常数，分别取0.09和0.4。

表1  验证算例边界条件参数

Tab.1  Boundary condition parameters of validation example

位置	边界条件	设置情况
入口	velocity_inlet	u和k采用风洞试验数据拟合
出口	outflow
顶部和两侧	symmetry	自由滑移
底部	wall	标准壁面函数

图3为本文的模拟结果与风洞试验结果^［6］以及其他学者的RANS方法模拟结果^［11］的对比情况。图中平均风速剖面经过归一化处理，纵轴为z/h，h为粗糙块高度 0.15 m。横轴为$u_{\mathrm{i}\mathrm{n}\mathrm{l}\mathrm{e}\mathrm{t}}$/u_ref，u_ref设置为$\mathrm{3}\mathrm{m}·{\mathrm{s}}^{-\mathrm{1}}$。从图3中可以看出，本文的计算结果和风洞试验结果以及其他学者的模拟结果非常接近。观察建筑平均高度H以下范围，与风洞试验结果拟合较好，说明质量较好的四面体网格能够还原实际的风场情况，并且平均风速剖面的模拟结果较好。

图3  平均风速剖面对比

Fig.3  Comparison of mean wind speed profile

2 不同建筑环境下风场的模拟

2.1　模拟工况设置

陈泂翔^［7］通过资料归纳和数据统计，指出我国大型城市的建筑平均高度变化范围为25~45 m，建筑密度λ（建筑平面面积与占地面积之比）的变化范围为0.15~0.35，建筑迎风面指数$f$ （建筑迎风面面积与建筑占地面积之比）的变化范围为0.15~0.45，并介绍了采用建筑平均指标来衡量整个城市建筑群特征的形态学方法。杨淳^［3］通过分析试验数据指出当来流经过超过8 km的城市建筑群，风剖面才会稳定下来。因此本文采取沿流向长度超过8 km的城市建筑模型，对具有大型城市建筑环境参数的模型进行系统性模拟。

对图4给出的模型参数，本文将建筑模型底边长B设置为25 m，建筑平均高度H设置为6~35 m，建筑密度$\lambda$设置为0.06~0.51，建筑迎风面指数$f$ 设置为0.03~0.71。

图4  建筑模型参数设置

Fig.4  Parameters of building model

根据图4中各模型参数可以计算出各个工况的建筑密度$\lambda$和建筑迎风面指数$f$。

$$\lambda =\frac{B\times B}{\left(B+D\right)\times \left(B+D\right)}$$

（1）

$$f=\frac{B\times H}{\left(B+D\right)\times \left(B+D\right)}$$

（2）

从公式（1）和（2）容易看出，固定底边长B后建筑迎风面指数和建筑密度以及建筑平均高度是相关联的，因此在工况设置时，主要考虑建筑密度和建筑平均高度的变化。

相比于陈泂翔^［7］提出的大型城市建筑环境参数的取值范围，本文还考虑了建筑高度不高，但是粗糙地貌长度较长的情况。本节对粗糙地貌最后一排粗糙元上方50~350 m高空处的风剖面采用指数律进行拟合。由于模拟工况较多，相同大小的计算域不能适用于所有工况。本文通过建筑间距D来控制计算域的大小。所有工况计算域的建立都遵循以下3个原则保证计算结果具有可比性：

（1）粗糙地貌的总长度必须超过8 000 m。

（2）建筑模型在垂直流向全部布置为5列，最外两列建筑距离计算域侧边界D/2长度。

（3）研究地貌距入口距离为500 m，距出口距离为2 000 m，计算域高度为1 000 m。该设置阻塞率小于3%，符合要求。

为了适应各种不同建筑模型，采用四面体非结构网格进行划分。图5为网格局部划分情况，所有工况网格总量在400万到1 000万之间。壁面设置一层1 m高度边界层网格，建筑附近的网格大小设置为5 m，满足城市风场模拟指南^［8］的要求。

图5  局部网格划分

Fig.5  Scheme of partial grid discrete

表2为所有工况基本参数设置。湍流模型采用Realizable k-$\epsilon$模型。求解设置压力速度耦合方式为SIMPLEC，动量方程和湍流模型方程非线性对流项采用second order upwind格式离散，压力插值格式采用second order，梯度插值方法采用least squares cell based。各变量和连续性方程的残差收敛标准设置为10^-8。由于建筑模型复杂，并且计算域较大，计算过程中出现了振荡收敛的现象，所以在每个工况计算过程中，监测目标风剖面10 m和300 m高度处风速变化情况，计算结束之前所有工况监测点的风速都已收敛。

表2  边界条件参数设置

Tab.2  Parameters of boundary condition

位置	边界条件	设置情况
入口	velocity_inlet	采用GB 50009-2012中的A类地貌平均风速剖面和AIJ 2004的湍动能和湍流耗散率
出口	outflow
顶部和两侧	symmetry	自由滑移
底部和粗糙元	wall	标准壁面函数

2.2　模拟结果分析

表3给出了所有工况的高度在50~350 m的平均风速剖面指数α的模拟结果，采用相关系数R的平方表示拟合效果，其中最小R²=0.976，表明拟合质量较好。图6为α随着建筑平均高度变化的规律。保持建筑密度不变，可以看出α随着建筑平均高度的增大而增大，并且增大的速度较快，基本呈线性增长，跨度达到了0.12左右。对于给出的4种建筑密度，可以看出建筑密度为0.11和0.17的工况α增长速度要大于建筑密度0.25和0.31的工况。

表3  平均风速剖面指数模拟结果

Tab.3  Simulated results of mean wind speed profile exponent

建筑平均高度/m	平均风速剖面指数
	$\lambda =$ $0.06$	$\lambda =$ $0.11$	$\lambda =$ $0.13$	$\lambda =$ $0.15$	$\lambda =$ $0.17$	$\lambda =$ $0.21$	$\lambda =$ $0.25$	$\lambda =$ $0.31$	$\lambda =$ $0.39$	$\lambda =$ $0.51$
6		0.25			0.27		0.28	0.27
9		0.27			0.29		0.29	0.28
15	0.28	0.30	0.30	0.30	0.30	0.31	0.30	0.30	0.29	0.27
18		0.31			0.31		0.31	0.31
25	0.32	0.34	0.33	0.33	0.35	0.34	0.33	0.32	0.32	0.30
30		0.35			0.36		0.35	0.33
35	0.37	0.38	0.37	0.38	0.38	0.38	0.37	0.36	0.35	0.34

图6  平均风速剖面指数随建筑平均高度变化规律

Fig.6  Relationship between mean wind speed profile exponent and average building height

由于我国建筑荷载规范（GB50009—2012）中对于平均风速剖面指数的规定，仅针对建筑平均高度进行判断，所以在图中将其与本文模拟的结果进行了对比。在各个平均高度和建筑密度的工况下，本文计算出的平均风速剖面指数α都要大于规范所给出的α值。规范整体α随着建筑平均高度的增长上升更快，但是由于是按照梯度增长，所以本文给出的分类方式要更加详细。比如建筑平均高度超出18 m的城市地貌，在我国很多大型城市中都有出现。而建筑平均高度较低，但是地貌范围较大的情况在城市中也比较多见，这些情况在规范中都没有涉及。

图7为建筑密度对平均风速剖面指数α的影响规律。当建筑平均高度保持不变，改变建筑密度时，平均风速剖面指数α也出现了一定的规律性。从图7中可以看出，平均风速剖面指数α随着建筑密度的增大，先增大再减小，这个现象在Macdonald等^［14］提出的对数律理论中也有提及。当建筑密度非常小，地面平坦，平均风速剖面指数自然很小。当建筑密度非常大，气流流入城市街区的部分较少，而大部分气流从建筑顶端掠过，此时风剖面相当于被抬升了一个高度，而风剖面的形状改变并不大，这很好地解释了平均风速剖面指数α的减小。说明真实的城市地貌中，对于建筑高度相近的城市街区，平均风速剖面指数随着建筑密度增大并不是一直增大，而是先增大后减小，并且使得α最大的建筑密度的区域落在0.1~0.3之间，这个峰值区域正是城市中比较常见的建筑密度范围。当然，如果建筑高度参差不齐，并且差距较大，可能需要另外讨论。

图7  平均风速剖面指数随建筑密度变化规律

Fig.7  Relationship between mean wind speed profile exponent and building density

2.3　建筑环境参数和平均风速剖面指数的关系

我国建筑荷载规范对于城市地貌类别的判断考虑的因素较少，并且划分类别不够细致，对α的数值估计偏小。本节采用模拟的城市风剖面数据，综合考虑了建筑平均高度和建筑密度两个参数对α的影响，给出了α关于这两个参数的公式。

前文已经分析了α随着建筑平均高度和建筑密度的变化规律，分别对两个参数的规律进行拟合。

$$f\left(x\right)=Ax+C$$

（3）

$$f\left(x\right)=a{\mathrm{e}}^{-\frac{{(\mathrm{l}\mathrm{n}x-b)}^{\mathrm{2}}}{c}}$$

（4）

式中：A、C、a、b、c为常数。

建筑平均高度的线性规律直接采用一次函数（式3）进行拟合。建筑密度的规律主要考虑两个特征：先升高后降低；峰值左右不对称。考虑采用对数高斯函数（式4）对其进行拟合。值得一提的是，建筑密度λ变化范围为0~1，当λ趋近于0或者1时，代表两种没有粗糙度的理想情况：地表为平面以及抬升了一定高度的平面，这两种情况的α都趋近于0。因此用对数高斯函数表达α的意义非常合理。

由于建筑平均高度相对于建筑密度对平均风速剖面指数α大小的影响更显著，因此在同时包含两个参数的拟合时，将建筑平均高度看作主要影响因素。于是得到平均风速剖面指数α关于建筑密度以及建筑平均高度公式的基本形式。

$$\alpha =\left(A\lambda +C\right)H+a{\mathrm{e}}^{-\frac{{(\mathrm{l}\mathrm{n}\lambda -b)}^{\mathrm{2}}}{c}}$$

（5）

采用梯度下降法对公式（5）中的参数关系进行回归分析，即找到一组常参数的解使得公式计算结果和实际模拟结果的均方误差$E_{\mathrm{N}\mathrm{M}\mathrm{S}}$最小。

$$E_{\mathrm{N}\mathrm{M}\mathrm{S}}=\frac{{\displaystyle \sum _{i=\mathrm{1}}^n} {\left(O_i-P_i\right)}^{\mathrm{2}}}{{\displaystyle \sum _{i=\mathrm{1}}^n} O_i{P}_i}$$

（6）

式中：O为实际模拟结果中的α值；P为采用拟合公式计算的α值。经过多次迭代计算，$E_{\mathrm{N}\mathrm{M}\mathrm{S}}$达到了0.000 539 38，此时相关系数R²为0.976。回归分析结果见表4。

表4  梯度下降法计算结果

Tab.4  Results of gradient descent method

参数	数值
A	-0.005 317
C	0.004 649
a	0.257 287
b	-1.200 340
c	10.727 3

工程中为了方便应用，将各参数保留两位有效数字。最终本文给出的平均风速剖面指数α关于建筑平均高度H和建筑密度λ的拟合公式为

$$\alpha =\left(-\mathrm{0.005}\mathrm{3}\mathrm{\lambda }+\mathrm{0.004}\mathrm{6}\right)H+\mathrm{0.26}{\mathrm{e}}^{-\frac{{(\mathrm{l}\mathrm{n}\lambda +\mathrm{1.2})}^{\mathrm{2}}}{\mathrm{11}}}$$

（7）

图8为α的公式预测值与模拟结果的对比，图中的数据点都集中在y=x的直线附近，可以看出公式（7）预测准确性较好。

图8  公式预测的平均风速剖面指数与模拟结果对比

Fig.8  Comparison of predicted mean wind speed profile exponent and simulated results

在反映建筑平均高度和建筑密度对平均风速剖面指数影响规律的图6和图7中，加入了采用公式预测的α值，得到图9和10。从中可以看出采用公式（7）预测的α值也能较好地反映建筑平均高度和建筑密度两个变量对α的影响规律。

图9  公式预测的平均风速剖面指数与建筑平均高度的关系

Fig.9  Relationship between predicted mean wind speed profile exponent and building mean height

图10  公式预测的平均风速剖面指数与建筑密度的关系

Fig.10  Relationship between predicted mean wind speed profile exponent and building density

3 结论

本文对具有不同城市建筑环境参数的模型进行了系统性的数值模拟，得出了充分发展的城市地貌的平均风速剖面指数模拟结果。通过对结果进行分析，得出了以下结论：

（1）其他条件固定不变，仅改变建筑平均高度H时，平均风速剖面指数α随着建筑平均高度H的增长呈线性增长。将其规律与规范进行对比，发现本文模拟得出的α值在各个建筑平均高度下都大于规范给出的α。α的模拟结果随建筑平均高度H的增长率要小于规范值，但是规范在建筑平均高度超过18 m后缺少更多的分类。

（2）其他条件固定不变，仅仅改变建筑密度λ时，平均风速剖面指数α随着建筑密度的增大先增大后减小。α达到最大时建筑密度的范围处于0.1~0.3之间，这个峰值区域正是城市中比较常见的建筑密度范围。

（3）从平均风速剖面指数α的变化区间来看，建筑平均高度的影响要大于建筑密度。综合考虑两个参数的影响，结合回归分析方法，给出了通过建筑环境参数计算平均风速剖面指数的经验公式。在工程应用中，对拟建建筑来流方向8 km的地块进行统计分析，得到建筑平均高度和建筑密度。将二者代入公式，可以得出平均风速剖面指数α，进而确定风速高度变化系数。

作者贡献声明

全 涌：确定研究方向，提供研究经费和专业指导。

曾一凡：进行数值模拟工作，对模拟得出的数据进行分析，撰写论文。

全 茜：图表及文字处理。