基于改进大涡模拟模型的风荷载验证

娄伶伶, 吴玖荣

(广州大学 广州大学-淡江大学 工程结构灾害与控制联合研究中心， 广东 广州　510006)



基于改进大涡模拟模型的风荷载验证

娄伶伶, 吴玖荣*

(广州大学 广州大学-淡江大学 工程结构灾害与控制联合研究中心， 广东 广州510006)

为了获得与近地大气边界层实际脉动风场特征相一致的三维脉动风速，文章首先应用DSRFG (Discretizing and Synthesizing of Random Flow Generation)方法，模拟生成计算流体动力学中大涡模拟算法所需的模拟入口边界脉动风速时程，然后采用改进的动态一方程模型作为大涡模拟的亚格子模型，计算长宽高为1∶1∶3的高层单体矩形建筑在0度风向角下的风荷载分布及特性，并与日本TPU风洞试验进行了对比.结果表明：①当高层建筑长宽比(B/D)=1时，扭矩的能量主要来自涡旋脱落而不是入口湍流的影响；②DSRFG结合改进的动态一方程模型能够较为准确地模拟得到高层建筑表面的平均及脉动风荷载特征.此方法有望应用到实际高层建筑风荷载评估中.

高层建筑; DSRFG; 大涡模拟; 亚格子模型; 风荷载; 风洞试验

高层建筑具有高、细、柔的特征，其自振频率较低且阻尼较小，过大的风速会使得高层建筑外围结构甚至主体结构损伤.由于能够获得风场作用下风荷载的脉动时程及相关统计信息，大涡模拟(Large Eddy Simulation, LES)已经成为计算流体动力学(Computational Fluid Dynamics, CFD)对高层建筑风荷载数值模拟研究的主要手段之一.能否较好实现LES在风工程中的应用，TAMURA[1]总结了3个关键点：①入口湍流风速的准确模拟；②亚格子湍流模型的选用及处理；③计算流体动力学方法中相关数值离散化算法的处理.

近地边界层风为具有高雷诺数的非定常流，且表现出较强的不均匀性和各向异性.因此，为获得与近地边界层实际脉动风场特征相一致的三维脉动风速，需要较为准确地模拟生成满足相空间相关性和湍流特性的入口边界脉动风速时程.目前CFD中主要采用2种入口边界脉动风速生成方法：①在辅助的流体计算域内通过recycling方法生成湍流风；②通过人工数值模型生成入口边界风.NOZAWA等[2]用第一种方法生成了入口边界.DSRFG[3]属于第二种方法，DSRFG方法可以生成满足特定频谱和空间相关性的脉动湍流场，且生成的脉动速度满足连续和自保持条件.

大涡模拟算法中采用的亚格子模型一般有如下几种方式：代数涡粘模型[4-5]、一方程模型[6]、两方程模型[7].然而上述亚格子模型在计算风工程实际应用中存在2个有待改进之处：①大多数亚格子模型适用于高精度数值格式和结构网格[8]，但是实际工程中结构比较复杂，一般采用精度相对较低的数值分析算法及无结构网格；②在高雷诺数流和网格相对粗糙情况时，代数涡粘模型的局部平衡假定与其表现出来的不平衡特征相矛盾.因此，在结合上述亚格子模型各自优缺点的基础上，有必要建立一种与实际工程CFD风荷载模拟需求情况较为接近的改进亚格子模型.HUANG[9]等采用结合动态一方程模型[6]和WALE模型[5]各自优点的改进亚格子模型，进行大涡模拟高层建筑风荷载，其特点有：①可用于高雷诺数流和粗糙网格的计算；②确定动态参数时不需要试验过滤，可用于低阶数值离散和无结构化网格；③用不同的机制处理亚格子能量的生成及GS(grid-scale)能量耗散.

本文采用DSRFG方法，模拟与高层建筑所处的近地大气边界层风场特征较为一致的入口边界风速时程，并采用改进的动态一方程模型[9]相对应的大涡模拟算法，计算了长宽高比为1∶1∶3的高层单体矩形建筑在0度风向角下的风荷载分布及特性，并与日本TPU风洞试验[10]结果进行了对比研究，分析LES获得的平均及脉动风荷载特性，验证了上述方法的有效性和准确性.

1　本文采用的大涡模拟相关算法概述

本文的大涡模拟CFD分析是在含有48个CPU的多处理器工作站进行的，该工作站可用于大规模的并行计算任务，采用计算软件为Flunet 6.3.

1.1计算模型与网格划分

采用的高层建筑模型[10]尺寸为40 m×40 m×120 m，几何缩尺比1∶400.CFD计算域尺寸：X、Y、Z3个方向为2 600 m×1 000 m×1 000 m.CFD计算时缩小模型风洞阻塞比不应超过3%[11]，本文模型的计算阻塞比为3%，满足模拟要求.

CFD计算时采用六面体结构化网格，图1为计算域网格的划分情况.为保证模拟精度并节省计算时间，对建筑物表面和模型附近的网格进行加密(图2)，并向外逐渐增大网格密度，划分出约70万个网格单元.

图1　CFD计算区域示意图

图2　模拟建筑物周边网格局部放大示意图

1.2CFD中入口边界脉动风速时程模拟

风工程中大多数实测的近地大气边界层脉动风速谱都服从冯·卡门(von Karman)谱[12]，为此本文采用DSRFG方法[3,13]，以模拟生成满足目标谱(冯·卡门谱)的CFD入口边界风速时程.其主要步骤如下：

(1)对于具有任意形式的3D脉动风速能量谱E(k)，均可由一系列离散的E(km)构成：

(1)

对于任一方向(此处假设为x方向)具有离散风速能量谱的脉动风速，其脉动速度场可由下式模拟：

(2)

(2)脉动风速在频域的合成：

(3)

(3.1)

(3.2)

ωm,n∈N(0,2πfm)，fm=kmUavg

(3.3)

其中，ζ，ξ∈N(0，1)，fm为频率，Uavg为平均风速，Ls为湍流积分尺度.a满足(0～1)之间的均匀分布.kmax=500，N=50时模拟结果与目标谱较为接近[12].因此，本文中kmax=500，N=50.

本文采用B类地貌，α=0.15.入口边界处平均风速度剖面采用如下的指数率形式：

(4)

其中，Z0=10 m，U0代表10 m高度处的平均风速，本文取15 m·s-1.

根据日本《AIJ建筑荷载建议》[14]，湍流强度表达式如下：

(5)

其中，zb为近地高度尺度，取7.5 m，zG代表梯度风高度取350 m.图3、图4分别为入口处平均风速、X方向湍流强度随高度变化示意图.

图3　平均风剖面

图4　X向湍流强度

图5给出了DSRFG法模拟的风场中的脉动风速功率谱与冯卡门谱的对比图，可看出2者吻合较好，说明入口湍流脉动风可以较为真实地模拟实际大气湍流边界层情况[3,13].

图5　DSRFG法生成的脉动风速功率谱与冯·卡门谱对比

1.3大涡模拟中改进亚格子模型概述

基于KAJISHIMA[6]的一方程模型，ksgs可以从方程(6)求得[9]:

(6)

(7)

最后一项εω是额外的耗散项，可以解释近壁面低雷诺数的影响.

根冠比是指植物地下部分与地上部分的鲜重或干重的比值，它的大小能反映出植株地上、地下部生长发育情况以及土壤的营养供应状况，该指标高则根系机能活性强，低则弱［22］。锦紫苏在不同栽培基质下，除处理⑥外，其他处理根冠比均高于对照组，说明不同配比基质的城市堆肥污泥有利于植物根系发达，增加植物的根系机能活性，从而有利于培养壮株，有利于作物的地上部分生长。试验证明，城市污泥堆肥对于锦紫苏生长具有好的适应性，并为污泥堆肥直接用于观叶类草本植物的种植提供了理论基础，在一定程度上促进了城市污泥堆肥的资源化利用。

(8)

本文所采用的改进亚格子模型具有如下优点：不需要试验过滤，可用于无结构网格计算，可用于高雷诺数流的模拟，采用Fluent UDF编程加以完成.

1.4CFD求解设置及边界设置

CFD中有限体积法对应的离散方程组的求解，一般采用收敛性较好的SIMPLEC算法[15].对一定偏斜的网格，Skewness Correction可减小收敛困难，默认为0.LES计算中，动量方程对流项的离散采用Bounded Central Differencing算法.时间项的离散采用二阶隐式方案.各项松弛因子均采用默认值.计算时间步长0.001 s，CFD模拟计算采用10 000步.

流场出口采用完全发展的出口边界条件，计算域顶部和两侧采用自由滑移壁面，计算域底面和高层建筑表面均采用不可滑移的壁面.

2　CFD数值模拟结果

日本TPU风洞试验包含了具有多种不同长宽高比值的矩形单体建筑模型，位于各类常用的地貌条件下的多测点高频测压风洞试验结果，为高层建筑结构抗风设计提供了在多风向角条件下的结构物表面详细的风荷载时程风洞试验数据[10].本文CFD模拟结果将与TPU数据库中长宽高之比为1∶1∶3的矩形单体高层建筑模型，在0度风向角下的风洞试验数据结果进行对比分析，以分析本文所采用大涡模拟CFD计算算法的有效性.

2.1平均风压系数

以本文所分析的高层建筑模型在Z=0.25 m、Z=0.15 m、Z=0.05 m 3个高度的部分测点为研究对象，各层布置的测点数与风洞试验相同[10]，为20个，测点布置见图6.

图6　测点布置图

(9)

图7为大涡模拟与风洞试验在模型3个高度处不同测点平均风压系数的对比图.1～5号测点在高层建筑迎风面，6～10号测点在右侧面，11～15号测点在背风面，16～20号测点在左侧面.

图7风洞试验与数值模拟各测点的平均风压系数对比

Fig.7Comparison of mean wind pressure coefficient by CFD and wind tunnel

通过对比分析可知：①在模型Z=0.25 m高度，迎风面为正压，呈现出中间大、边缘小(最大0.89，最小0.49)的变化趋势.5号测点与6号测点之间风压系数迅速降低，说明来流遇到建筑物的阻挡在模型两侧分离较强，边缘处负压梯度较大.侧风面与背风面均为负压，且各测点所得的模拟结果精度比较高，特别是10号测点至17号测点，误差均在10%以内.迎风面误差较大，可能是迎风面与湍流入口距离不足或边界条件设置的影响所导致，但是整体趋势吻合较好；②当来流平行于结构对称轴时，CFD数值模拟的平均风压基本为对称分布，而风洞试验数据对称性较差；③Z=0.15 m、Z=0.05 m时与Z=0.25 m时平均风压系数的分布类似.

2.2脉动风压均方根系数

脉动风压系数反映脉动能量的大小是脉动风荷载的重要特征.脉动风压均方根系数如下：

(10)

图8可见：①脉动风压均方根系数与平均风压系数具有相似的分布趋势；②来流在迎风面边缘分离再附着，导致迎风面脉动风压值较大，变化梯度也相对较大，在其他部位变化相对平缓；③随着高度的变化，建筑背风面脉动风压均方根系数差值在10%以内，说明背风面涡旋的作用比较均匀；④脉动风压系数随着高度的降低而逐步变小，说明建筑低处含有的脉动能量相对较小.

图8　数值模拟各测点的脉动风压均方根系数

Fig.8Root mean square coefficient of fluctuating wind pressure by CFD

2.3层脉动风荷载功率谱

高层建筑风荷载分为2个主轴方向上的顺风向、横风向荷载以及绕Z轴的扭矩.沿高度将建筑分为6层，各测点与风洞试验测点布置相对应.将CFD计算时各层测点每个时间步上的脉动风压，进行积分可得到各层的脉动风荷载时程数据，然后从频域的角度表达层风荷载的特性.限于篇幅，本文只取第2、4、6层进行对比分析.

2.3.1层顺风向力功率谱

①图9a、b、c分别为第2、4、6层的顺风荷载归一化功率谱，低频段(无量纲频率0.04～0.40)能量分布较均匀，随后谱密度开始下降，但是衰减速度比风洞试验结果慢，可能是因为CFD模拟中入口处与模型迎风面的距离比风洞试验小；②能量主要集中在低频段，且与图5所示的入口脉动风速功率谱形状相似，说明入口风湍流是产生顺风向风荷载的主要原因；③图9d为前述3层数值模拟功率谱对比，在低频段第6层功率谱值最小，第2层最大，因此，低频段较高层的功率谱能量比低层小，可能是脉动风在建筑迎风面下沉，在地面形成回流，湍流较大导致的.高频段不同层上没有明显的变化.

2.3.2层横风向力功率谱

①图10a、b、c分别为第2、4、6层的横风向荷载归一化功率谱，低频段(无量纲频率0.04～0.10)能量高于风洞试验结果.高频段与风洞试验结果较为接近；②不同层功率谱曲线均在无量纲频率0.1处出现一个明显的尖峰，其主要是涡旋脱离造成的；③由图10d知，低频段第2层峰值最大，第6层峰值最小，高频段3者数值基本接近.

2.3.2层扭转风荷载功率谱

①图11a、b、c分别为第2、4、6层的扭转风荷载归一化功率谱，当频率小于无量纲频率1时模拟值结果大于风洞试验值，大于1时模拟值小于风洞实验值，在无量纲频率0.10附近，出现一个明显的尖峰，与横风向功率谱峰值对应的频率接近，涡旋脱落现象明显；②扭转谱与风洞试验结果的误差，与层顺风向及层横风向谱和相应风洞试验结果对比来看，说明本文CFD模拟结果对层扭转谱的模拟精度要相对低一些.但CFD数值模拟结果所反映的扭转谱的涡旋脱落现象，与风洞试验结果基本相似.

图9　层顺风向荷载功率谱

图10　层横风向荷载功率谱

图11　层扭转风荷载功率谱

3　结　论

本文应用DSRFG方法，模拟生成CFD大涡模拟算法所需的入口边界脉动风速时程，采用改进的动态一方程模型作为大涡模拟的亚格子模型，计算了长宽高为1∶1∶3的高层单体矩形建筑在0度风向角下的风荷载分布及特性，并与日本TPU风洞试验进行了对比分析，得到以下结论：

(1)DSRFG可以生成大涡模拟算法中与近地大气边界层实际脉动风场特征相一致的入口边界脉动风速场.DSRFG方法结合改进的动态一方程亚格子模型，运用于CFD大涡模拟算法中，可以较好地预测高层建筑表面的风荷载，数值模拟与风洞试验数据结果基本一致.

(2)LES模拟得到的高层建筑表面平均风压系数与风洞试验数据吻合较好，反映了典型钝体绕流的特征；来流平行于结构对称轴时，风洞试验的平均风压对称性很差，而数值模拟的对称性很好，体现了数值模拟的优势.

(3)LES的脉动风压均方根系数与平均风压系数具有相似的分布趋势.随着高度的增加，迎风面及两侧面脉动能量有逐步增大的趋势，背风面涡旋作用较均匀，没有明显变化.

(4)LES对层风荷载归一化功率谱的模拟表明，入口湍流风是产生顺风向荷载的主要原因；涡旋脱落是造成横风向荷载的主要原因，功率谱曲线在无量纲频率0.1处出现明显的尖峰；扭矩一般由入口风湍流与涡旋脱落相互作用造成，但对于本文的模型B/D=1，扭矩主要来自涡旋脱落的影响.

[1]TAMURA T. Towards practical use of LES in wind engineering[J]. J Wind Eng Indust Aerod, 2008, 96:1451-1471.

[2]NOZAWA K, TAMURA T. Large eddy simulation of wind flows over large roughness elements[C]∥Proceedings of the 4th European and African Conference on Wind Engineering (EACWE4),2005, Prague, Cezch Republic.

[3]HUANG S H, LI Q S, WU J R. A general inflow turbulence generator for large eddy simulation[J]. J Wind Eng Indust Aerod, 2010, 98(10/11):600-617.

[4]GERMANO M, PIOMELLI U, MOIN P, et al. A dynamic subgrid-scale eddy viscosity model[J]. Phys Fluid, 1991, 3(7):1760-1765.

[5]NICOUD F, DUCROS F. Subgrid-scale stress modelling based on the square of the velocity gradient tensor[J]. Flow Turbul Combust, 1999, 62(3):183-200.

[6]KAJISHIMA T, NOMACHI T, KAJISHIMA T, et al. One-equation subgrid scale model using dynamic procedure for the energy production[J]. J Appl Mech, 2005, 73(3):368-373.

[7]GALLERANO F, PASERO E, CANNATA G. A dynamic two-equation sub-grid scale model[J]. Contin Mech Thermod, 2005, 17(2): 101-123.

[8]TSUBOKURA M. Proper representation of the subgrid-scale eddy viscosity for the dynamic procedure in large eddy simulation using finite difference method[J]. Phys Fluid, 2001, 13(2):500-504.

[9]HUANG S, LI Q S. A new dynamic one-equation subgrid-scale model for large eddy simulations[J]. Intern J Numer Method Eng, 2009, 81(7):835-865.

[10]林宁，梁波，田村幸雄. 高层建筑层风力特性实验研究[J]. 振动工程学报，2003(12): 409-414.

LIN N, LIANG B, YUKIO T. Experimental investigation on local wind force of high-r ise building[J]. J Vibr Eng,2003(12):409-414.

[11]黄本才. 结构抗风分析原理及应用[M].上海：同济大学出版社, 2008.

HUANG B C. Principle and application of structural wind resistance analysis[M]. Shanghai: Tongji University Press, 2008.

[12]LI Q S, WU J R, LIANG S G, et al. Full-scale measurements and numerical evaluation of wind-induced vibration of a 63-story reinforced concrete tall building[J]. Eng Struct, 2004, 26(12):1779-1794.

[13]熊涛, 吴玖荣. 大涡模拟中入流边界脉动风速模拟相关参数取值分析[J]. 广州大学学报：自然科学版, 2014, 13(3):60-67.

XIONG T, WU J R. Parameter analysis for modeling the inlet fluctuated wind velocity in large eddy simulation[J]. J Guangzhou Univ:Nat Sci Edi,2014,13(3):60-67.

[14]Architectural Institute of Japan. Recommendations for loads on buildings[M]. (English Version). Tokyo: AIJ Press, 2004.

[15]VAN DOORMAAL J P, RAITHBY G D. Enhancements of the simple method for predicting incompressible fluid flows[J]. Numer Heat Transf Appl, 1984, 7(2):147-163.

【责任编辑： 周全】

Verification and evaluation on wind loading of large eddy simulation by the improved sub-grid scale model

LOU Ling-ling, WU Jiu-rong

(Guangzhou University-Tamkang University Joint Research Center for Engineering Structure Disaster Prevention and Control, Guangzhou University, Guangzhou 510006, China)

In order to obtain the time history of 3-dimensional fluctuating wind speed with the same statistical characteristics as the wind of near ground atmosphere boundary layer, an inflow turbulence generation method called DSRFG (discretizing and synthesizing of random flow generation) is applied to generate the time record of inflow turbulence wind speed for Large Eddy Simulation in CFD, meanwhile an improved dynamic one-equation model is adopted as the sub-grid scale model for Large Eddy Simulation method. These two methods are combined for the CFD simulation on the wind load distribution and characteristics of an isolated rectangular tall building with 1∶1∶3 ratio (the ratio between length, width and height) in zero degree wind direction. The CFD numerical simulation results are further compared with those experimental wind tunnel test results on the same building from TPU wind tunnel laboratory in Japan. All these comparisons indicate that: the fluctuating kinetic energy of wind-induced torsional moment mainly results from the influence of vortex shedding other than the inlet turbulence for the test building; the combination of DSRFG with the improved dynamic one-equation model in the large eddy simulation method in CFD can accurately obtain the value and their characteristics of mean and fluctuating wind loads acting on the isolated rectangular tall building, which jointly testify the feasibility of its application to the design wind load evaluation for actual high rise buildings.

tall building; DSRFG; large eddy simulation; sub-grid scale model; wind load; wind tunnel test

2016-02-22;

2016-04-12

国家自然科学基金资助项目(51378134，51578169)；广东省高层次人才资助项目

娄伶伶(1989-)，女，硕士研究生.E-mail：loulinghcool@163.com

. E-mail: jrwuce@qq.com

1671- 4229(2016)04-0050-07

TU 312.1

A