airpak-建筑群流场模拟中边界条件的设定

建筑群流场模拟中边界条件的设定
摘要：CFD(Computational Fluid Dynamic)模拟技术近几年来随着计算机性能的大幅提高和模拟技术的进步，越来越得到广泛重视和应用。

边界条件的设定是直接关系到模拟成败的关键。

本文对边界设定中几个常见的问题分别进行了分析。

计算区域根据实际情况选取,水平方向参考取10倍最高建筑的高度，高度方向取3-4倍建筑物高度。

计算区域的顶面和侧面对于对称边界，速度边界和壁面这3种边界可以任意选择。

出流面可设为完全发展出流边界或速度边界。

对于大气层结中性情况，给出了风速廓线的计算公式，并介绍了一种为节约计算能力而分阶
段划分网格的方法。

关键词：边界条件建筑外流场数值模拟
1 引言
风是构成室外环境的重要因素。

建筑群造成的流场变化会对其周围环境会产生一定的影响，特别是高层建筑更容易引起建筑周围行人高度局部风速的增大，产生不良的风环境，造成人的不适，降低环境质量。

如果在建筑规划设计阶段进行模拟研究，问题就可以得到预防和解决。

另外流场模拟对于结构风荷载以及自然通风等领域的研究也有价值。

风洞模型试验和实测是以往建筑流场研究的主要手段，但都费时费力。

CFD(Computational Fluid Dynamic)模拟技术近几年来随着计算机性能的大幅提高和模拟技术的进步，越来越得到重视和广泛应用，并取得了不少的成果。

作为有力的工具，数值模拟也已被用于建筑内外的流场研究中。

在模拟的过程中，除了湍流模型的选择，网格的划分等因素，边界条件也是直接关系到模拟成败的关键环节。

目前，数值模拟技术还处于发展阶段，只有不断与风洞试验和实测的结果进行比较，总结经验，对各方面进行改进，才能使之更加成熟完善。

本文结合模拟的例子，对边界设定中几个常见的问题分别进行了分析讨论。

2 模拟例子的情况简介
本文通过结合一个模拟的例子来说明一些有关边界条件方面的问题。

例子是模拟上海某繁华商业街上一个供行人休息的小广场周围的行人高度（1.5至2米）的流场，结合日照等因素综合研究建筑对人及人的行为的影响。

小广场周围的大楼高度都在100米左右，形体及其造成的流场情况复杂。

建筑分布图如图1所示：
图1 广场周围建筑平面图
城市中下垫面结构复杂。

因在同一高度上增温冷却情况不均一，往往引起热力湍流，在高温区有上升气流，在低温区有下沉气流。

热力环流只在盛行风速极小或静风时才表现出来[文献1]。

因此，在本模拟中可忽略辐射及浮升力的影响。

由于建筑的存在而改变了流动结构，流动不能充分发展。

标准k-e模型对于某些情况过大估计了k值，但对于速度场，标准k-e模型与风洞试验结果较为接近[文献2，3]，此外其波动小，有较高的精度和效率，在低速湍流工程上使用广泛。

本文选用标准k-e湍流模型，近壁面则采用标准壁面函数法计算。

3 梯度风的边界模拟
任一时刻风速都是随高度变化的。

在摩擦层中，对于风随高度的变化规律，主要考虑摩擦力的影响。

风吹过地面时，受到地面上各种粗糙元(草、庄稼、森林、房屋建筑等等)产生的摩擦阻力作用而使风的能量减少即风速减小，这一层受地球表现摩擦阻力影响的大气层称为大气边界层，它的厚度约为1千米，建筑物模拟的计算区域一般不会超出这个高度。

对于垂直方向风速分布的描述常用以
下两类公式：
3.1 指数律
3.1.1 假定大气边界层风速剖面符合简单的幂指数分布规律，指数a在梯度
高度δ内保持不变，而δ本身只是a的函数，即：
(1)
式中Vo为在Zo高度的平均风速。

参考高度Zo一般采用10 米。

Vo即取
10米高度的平均风速。

当Z大于某一高度Z
H 时，风速趋于均匀V=V(Z
H
),Z
H
就是
地转风出现的高度。

不同的地面条件，幂指数a不同，一般来说，a在0.14-0.40内取值。

3.1.2 指数a的选取
对于指数a的值各研究者提供的值差异很大。

例如，波锐辛科等曾研究过在近地面100米高度以下，指数a与粗糙长度Zo（计算方法见下文公式4，5）具有线性关系：
a=0.12Zo+0.18 (2)
M.Hussain [文献4] 指出指数也可以用a=Ae-bv式得出，式中A,b是由当地实际条件取得的经验值。

结果与当地实测结果相符。

我国“建筑结构载荷规范”GBJ 9-87中不同类型地表面下的a值与梯度风高度（即大气速度边界层厚度）的关系如表1所示：
地面类型适用区域指
数a
梯度风
高度
A近海地区，湖岸，沙漠地区0.1
2
300m
B田野，丘陵及中小城市，大
城市郊区
0.1
6
350m
C有密集建筑的大城市区0.2400m
3.2 对数律
3.2.1 对数律是由理论推导得出的。

具体推导见[文献5] ：
(3)
式中表示在地转风高度以下，摩擦层内Z高度出的平均风速；V
*
表示摩擦
速度；k表示卡曼常数，一般取0.4；Z表示高度；Z
表示粗糙度长度（在中性稳定平衡下，垂直风速廓线外延到风速等于零的高度）。

上式的适用范围可由
0.02 V
*
/f确定，强风时在200米的范围内能满意地表达风速分布。

粗糙度长度Z
随地点的不同，相差很大，从零点零几到几米不等。

介绍两
种常用的经验公式：
1）（4）
式中V
1,V
2
分别为Z
1
,Z
2
高度的平均风速。

只要测出两高度及两个相应的风速，
就能算出粗糙度长度Z
0。

2）（5）
式中H表示障碍物高度（米）；a表示气流接触到的障碍物侧面积（米2）；
A表示障碍物覆盖面积。

常见的地形Z值如下表2所示：
海面海滨矮草地，干旷草
原田野休耕地
森
林
0.00 1
0.00
50.01-0.040.05
0.02-0.0
3
2.
松树林
果树
林稀疏建成市郊
密集建成
市郊
大都市中
心
0.9-
1.0 1.00.2-0.40.8-1.2
2.0-
3.0
3.2.2 对数律存在的问题
从公式3可以看出，对数律在近地面小于粗糙长度Z
时求得的风速为负值，这显然是不合理的。

对于研究建筑风荷载，这点问题尚可以忽略，但对于行人高
度的流场情况研究，这是无法忍受的。

纽陶就Ekman边界层内风速对数廓线形式提出以下公式:
（6）
式中表示在Z高度的平均风速；V
G 表示地转风速，V
G
=，式中为
空气密度；f为科氏力参数，可根据所在地的纬度查得；为该地的水平压力
梯度，从气象台可查得；c表示阻力系数，c阻力系数与无量纲的罗斯贝数（Rossy Number Ro）有关。

Csanady给出其关系式为：
经观测和模拟试验证实，上述公式符合实况。

与公式3相比，公式6避免了风速出现负值的情况。

笔者认为在有条件取得地转风的大小的情况下，可以用上
式代替对数律使用。

3.3 对数律与指数律的比较和共同存在的问题
目前，对于近地面风速的描述到底是用对数律[文献6]还是指数律[文献7]尚无定论，两种公式都有实测支持。

还有待微气象方面的学者进行进一步的观测结果来验证。

目前大部分使用CFD工具的学者倾向于使用指数律。

值得注意的是对于不同的大气层结情况，风速廓线是不同的。

只有在中性时，才近似符合对数律或指数律。

因此对于夜晚等大气层结非中性情况，使用指数律或对数律描述风速会有很大的偏差。

平均风速与对数廓线的偏差如表3所示：
情况1情况2情况3
高度502001515
不稳定层
结-1%-4%-4%≈-22%
稳定层结1%4%2%≈12%
注：情况1，2，3，分别为10米高度时均风速约25m/s,12m/s,5m/s时的情
况。

表3数据表明情况3偏差达10%以上。

需根据实际情况进行修正。

但是迄今为止，尚无可用于任意层结条件下的边界层风场的公式，虽有些工作计入了层结，但对层结处理过于简单且主要用于理论研究。

3.4 本例中的设置
本例根据上海市气象资料冬季主导风向为北风和西北风，平均风速为3m/s。

对于舒适性研究风速不应只用平均风速，这样的结果偏保险。

最常见的最大风速等于平均风速的1.1倍[文献8]。

采用指数律，考虑到地处市中心和前方的绿地，
指数a值折衷取0.35，风速廓线为：，风向西北。

4 计算区域边界条件的设定
4.1 顶面和侧面的边界设定
在无法对边界进行实测的情况下，可参考风洞模型试验的做法，将建筑群完全浸没在风里。

整个计算区域就如同一个大风洞，计算区域的顶面和侧面相当于风洞的壁面。

笔者在查阅一些文献后发现，顶面和侧面有3种边界可以选择：对称边界，速度边界和壁面。

这三种边界都有人使用，有将顶面和侧面都设成三者之一的，也有分别使用其中两种的。

其中速度边界设置为k,ε和速度三项的梯度为零，法向速度为零。

笔者分别将边界设成以上三种情况进行了模拟，三种情况1.5米高度截面速度矢量图如图2，图3，图4所示：?
图2 使用对称面边界图3 使用速度界面边界图4 使用壁面边界
从结果可以看出，除了靠外的一些部分流场有所不同，在内部主要关心的区域计算结果是基本相同的。

以上三种边界对内部流场的影响不大，即对最后结果影响不大。

虽然三者形式不同，共同的特点是通过对称平面没有扩散流量。

在建筑与计算区域各面之间有足够距离的情况下，三者是等效的，任取一种都可以。

4.2 入口面的边界设定
对入口边界的设定没有异议，笔者查阅的文献中都设为速度边界。

但是入口
湍流条件的设定意见并不统一。

4.2.1 描述大气湍流的最简单的参数是湍流强度。

[文献5]给出如下的湍流
强度剖面：
式中V
*表示摩擦速度；V
（z）
表示Z高度的平均速度；一般假
设β不随高度变化，对于不同粗糙长度的β值
Z
0.0050.070.31 2.5β 6.56 5.25 4.854
4.2.2 [文献9]根据上海和广州的实测提出如下的湍流强度剖面。

4.2.3 [文献10] 给出如下的湍流条件：
式中V
*表示摩擦速度；Z表示高度;常数C
μ
=0.09；
卡曼常数k=0.4。

笔者将几种条件带入，对计算结果进行比较，发现入口湍流条件，对结果影
响不大。

4.3 出流面的边界设定
对于出流面的设定，可设成完全发展出流边界，即各项的梯度都为0，这就要求在建筑群的下游有足够的长度；或设成与入口相同的速度边界，当然是基于风速在流过建筑后不变的前提，水平剪切使边界层内流动减速，因而边界层增厚，但在大气边界层内，在梯度高度之下，水平压力梯度只是部分地被科里奥利力平衡，它将再赋能量给流体，从而抵消边界层的增长，因此流动的水平均匀性就能够保持[文献5 ]。

笔者试验结果是两者对结果的影响不大，使用两种边界的速度等值线图如图5,图6所示。

可任选其一。

图5 使用充分发展边界的剖面图
图6 使用速度边界的剖面图
4.4 地面的边界设定
路面和建筑的墙面可以设成壁面。

然而除了建筑本身，周围栽种的植物尤其是枝叶茂盛的树木对近地面流场风速的阻碍作用很大，种植枝叶茂盛的行道树的人行道的风速比街道中心风速减低20-30%；公园浓荫下的风速为街道中心风速的50%[文献1]。

如不进行修正，风速结果可能偏大，对于行人高度的舒适性研究这点事不能忽略的。

由于植物形态的复杂性和多样性，对其进行建模是很困难的。

[文献11] 对树木对流场影响的模型进行研究，使用图7 所示薄片组成的树
木模型。

图7 树木模型
将抗力项F
i =-C
fi
<U
i
>加到运动方程的右边，将F
k
=<U
i
>F
i
加到k的
输送方程右边，将 F
ε= (L=) 或 F
ε
=加到ε的输送方程右
边，进行修正。

其中C
fi 是抵抗系数，a 是叶面积密度，C
pε
是模型系数。

另一种
可以考虑的简便方法是结合经验按照树木的高度以及繁茂程度等因素，对不考虑树木影响的情况计算得出的风速乘以一个适当的削减修正系数。

5 计算区域大小的选择及一种计算方法
5.1 计算区域大小的选择
计算区域指的是建筑外的立方体（当然也可根据实际建筑的特点用其他形状）。

它的大小将直接影响模拟的准确性和计算的时间。

由于计算机性能的限制，计算区域大小及与其相关的网格数不可能太大。

实际上，太大也没有必要，因为靠外的流场是我们不关心的，计算区域太大反而为计算造成困难。

根据有关文献和笔者的经验，计算区域通常长度方向取10倍建筑的高度，高度方向取3-4倍
建筑物高度已能满足要求。

不同高度长度深度建筑背后的涡旋区范围见[文献1]，虽然文中列出的只是单幢建筑物的情况，另外建筑群对周围流场的影响更加复杂，不但与高度长度深度有关，而且与周围建筑的相互位置以及建筑群迎风方位有关，但大致上可以认为前方影响3倍高度，后方影响6倍高度，侧面影响2倍高度，如周围有建筑物遮挡，影响的距离减小。

根据此原则，可将周围建筑分类选择，除去对关心区
域没有影响的建筑。

图8，图9，图10分别是本例计算区域高度为200米，300米与400米三种情况的垂直剖面风速等值线图。

结果看出200米时，流场的发展受到计算区域的阻碍，300和400米时就逐渐基本不受阻碍了，差别不大了。

考虑到计算机的性
能，本例中计算区域高度选用300米。

图8 高度200米风速等值线图
图9 高度300米风速等值线图
图10 高度400米风速等值线图
建筑群后留出的区域在计算机性能允许的条件下应尽量大，使得流场能够充分发展，与假设的边界条件相符，一般要有最高建筑的6倍。

建筑群前方留出的空间可以相对小一些，但因为在风吹到建筑物时，在滞止点部分风会向下吹，在建筑前方造成一个旋涡区，对建筑前的流场产生影响，所以不宜小于3倍建筑高度。

由于建筑对侧面流场影响范围较小，宽度方向的空间只要满足计算区域两边的面到各建筑有各建筑高度的2倍即可，使得侧面和顶面边界对计算区域的流场
不会产生太大的影响。

国外研究表明建筑模型在整个风洞剖面中控制在1-2%左右最为理想，最大不要超过3%[文献12] 。

由于数值模拟与风洞实验的原理相似，这个原则也可用在计算机模拟中来检验计算区域大小的选择。

但笔者实际感觉没必要达到以上的
要求。

当然，最好根据经验和实际情况（包括计算机的性能），进行几次试算来最终确定合适的计算区域。

一般先算只有最近的建筑的情况，然后由近到远添加一些可能产生影响的建筑，将前后两次计算的结果进行对比。

如果后一次的结果与
前一次相比变化不大，就意味着不必在这个方向添加建筑了。

在各个方向试算几
次后就能确定较理想的计算区域大小。

5.2 本例的计算区域选择
由于广场南面为淮海公园的大片绿地，入口前也有大片的绿地，本例中都作为空地。

本例的计算区域为1400mX1000mX300m。

如图11所示：
图11 计算区域
5.3 一种分阶段划分网格的计算方法
对于数值模拟，计算时间和计算能力是非常宝贵的，尤其在计算机性能不够的条件下。

下面介绍一种既节约计算时间又能较好的完成任务的方法文献[14,15,16] 。

过于粗糙的网格（受计算机性能限制，本例中,壁面附近网格间距4米，远离壁面的空间网格间距15米）使流场的一些特性无法表现，模拟的准确性受到影响。

由于大部分计算区域的流场是我们不关心的，所以可以先用粗糙的网格对整个区域进行计算，然后用细致的网格对关心的重点区域（本例中是广场周围150mX150mX50m
的区域，细网格间距1
米）进行重新划分。

最后用第一次计算结果中的速度值作为第二次计算的边界条件进行计算，以求得较为精确的流
场。

图12 粗糙网格下1.5米高度速度等值
线图图13 细致网格下1.5米高度速度等值
线图
从两张图的比较可以看出，建筑附近流场基本是相似的，后一张的结果更加细致，但考虑到前一张的结果使用的是非常粗糙的网格，进行第二次用细网格计算是完全有必要的。

其他截面的情况也基本和1.5米截面相同，由于篇幅有限就不贴出了。

虽然只读入速度场的信息，但由于速度起主要作用，所以这种做法是合理的。

从结果看，可以认为这种做法是有效可行的。

6 结论
1) 计算区域根据实际情况选取,水平方向参考取10倍最高建筑的高度，高
度方向取3-4倍建筑物高度。

2）只要空间足够大，计算区域的侧面和顶面3种边界：对称边界，速度边
界和壁面都可以选择，对结果没有影响。

3）入口界面设成速度界面，湍流条件对结果影响不大，湍流强度可取15%。

出流面可设为完全发展出流边界或速度边界。

4）风随高度的变化可参考文中给出的式子模拟。

5）进行两次网格划分和计算来节约计算时间的方法是合理并且可行的。

参考文献：
1.周淑贞，张超.城市气候学导论.上海：华东师范大学出版社，1985
2.王嘉松，黄震.上海交通大学学报，2002,36(10).1496-1499
3.杨伟,顾明.同济大学学报，2003.6:647-651
4.M.Hussain. Energy Converstion and Management,2002,43:467-472
5.埃米尔.希缪,罗伯特 H 斯坎伦. 风对结构的作用---风工程导论.上海:
同济大学出版社,1992
6.陆龙骅，卞林根，程彦杰.应用气象学报,2002,1:34-42
7.植石群，钱光明，罗金铃.热带地理，2001,6:131-134
8.B.B.巴图林.工业通风原理.北京：中国工业出版社，1965
9.田浦.第三届全国结构风效应学术会议论文集,1988:81-86
10.倫裕發, 持田灯, 吉野博.建築環境工学研究年報2001年度Vol
9:309-310
11.加藤敦子, 持田灯, 吉野博.建築環境工学研究年報2002年度Vol
10:301-302
12.关滨蓉,马国馨.建筑学报,1995,11:44-49
13.Qingyan Chen. Proceedings of the 2nd International Workshop on Energy and Environment of Residential Buildings in China:142-151
14.Shuzo Muracami, Akashi Mochida, Shisuke Kato. 建築環境工学研究
年報2002年度Vol 10:405-426
15.Yoshihide Tominaga,Akashi Mochida,Hiroshi Yoshino.建築環境工学
研究年報2002年度Vol 10:281-289
. .。