基于两相流模型的大跨屋面风吹雪数值模拟分析

Numerical Simulation Analysis of Snowdrift on Large-Span Roofs Based on a Two-Phase Flow Model

摘要

风吹雪引起的积雪不均匀分布是大跨屋面雪灾下破坏的主要原因之一，现有规范中未考虑风速、风向对屋面雪荷载分布的影响，本文基于两相流模型，对某大型体育场屋面风致雪漂移进行数值模拟，研究不同风速、风向下屋面积雪的不均匀分布情况，讨论屋面形状对积雪分布的影响。结果表明：风速不变时，位于屋面背风侧的雪颗粒发生沉积，侧风面的雪颗粒被侵蚀；低风速下，屋面积雪分布趋于均匀，受来流风向影响较小；高风速时，屋面积雪不均匀分布明显，在不同风向下，同一区域可能处于完全侵蚀或积雪沉积两种不同状态；随风速增大，屋面积雪侵蚀面积增大，相比日平均风速，50年重现期风速作用下屋面积雪总量减少1/3，但局部最不利雪压增大两倍以上。大跨屋面的雪荷载分布受风速、风向影响较大，结构设计时有必要综合考虑不同风速、风向对雪荷载的影响。

Abstract

Unbalanced snow loads caused by snowdrift are one of the main reasons for damage to large-span roofs under snow events. Existing codes do not consider the effect of wind speed and direction on snow load distribution on roofs. In this study， numerical simulations of snowdrift on a large stadium roof were conducted based on a two-phase flow model. The uneven distribution of roof snow under different wind speeds and directions was investigated， and the influence of roof shape was discussed. The results show that when the wind speed is constant， snow particles located on the leeward side of the roof tend to deposit， and snow particles on the sidewind side tend to erode. At low wind speeds， the distribution of snow on the roof tends to be uniform and less affected by wind direction. At high wind speeds， the uneven distribution of snow on the roof is obvious， and in different wind directions， the snow particles in the same area may be in two states： complete erosion or deposition. As wind speed increases， the area of snow erosion increases. Compared with the snowpack under the daily average wind speed， the snowpack under the 50-year return period wind speed is reduced by 1/3 and the most unfavorable snow pressure is more than twice as large. The snow load distribution of large-span roofs is greatly affected by wind speed and direction， and it is necessary to comprehensively consider their impact on snow loads during structural design.

Graphical abstract

关键词

风致雪漂移 / 大跨屋面 / 两相流模型 / 数值模拟 / 积雪不均匀分布 / 侵蚀 / 沉积

Key words

snowdrift / large-span roof / two-phase flow model / numerical simulation / uneven snow distribution / erosion / deposition

某大型体育场长约300 m，宽约250 m，最高约63 m，位于西安新中心新轴线核心位置^［1］。该建筑为大跨度空间结构，结构跨度大、外形不规则，且平立面布局较为复杂（图1），与普通建筑相比，此体量结构屋面更易受到风吹雪不均匀分布的影响。现行《建筑结构荷载规范》（GB 50009—2012）^［2］（简称：荷载规范）中没有明确规定复杂屋面的积雪分布系数，为获取合理可靠的雪荷载分布，有必要对风吹雪下屋面雪荷载分布进行模拟预测。

屋面的风吹雪数值模拟研究多将雪颗粒看作连续介质进行模拟，即采用Euler-Euler方法来模拟。BEYERS等^［3］基于Euler-Euler方法采用两相流模型模拟了立方体周边风致雪漂移，并与实测结果进行了对比，数值模拟结果和实测结果吻合较好，证明数值模拟可以用于风吹雪模拟。孙晓颖等^［4-5］以防雪栏和立方体周边雪分布为研究对象，采用Mixture多相流模型，对比了不同湍流模型对模拟结果的影响，并对鞍形屋面、球形屋面等4种典型大跨度屋面风致雪漂移进行了数值模拟，探讨了不同形式屋面的积雪分布规律。刘博雅等^［6］分别采用风洞试验和数值模拟对高低屋面积雪分布进行研究，验证了Realizable k-\epsilon湍流模型的可靠性，结果表明：与高屋面迎风时相比，低屋面迎风会导致屋面积雪分布更不均匀，且积雪分布系数可能超出荷载规范的建议取值。ZHANG等^［7］采用Realizable k-\epsilon湍流模型对阶梯屋面雪荷载进行数值模拟，得出阶梯屋面积雪分布系数最大可达3.44，低屋面最佳防雪坡度为11°。孙芳锦等^［8］采用Mixture多相流模型对4种大跨度双曲屋面在风雪流下的压力分布进行了模拟分析，得出风雪流下屋面的压力系数大于单独风作用下的压力系数的结论。大型场馆屋面几何造型复杂，无法归类于某种简单造型。周晅毅等^［9］、徐爱民等^［10］、卢春玲等^［11］分别对首都国际机场、长沙机场T3航站、吉林新火车站的复杂屋面进行了风吹雪模拟，计算分析了来流风作用下屋面积雪发生迁移后的分布情况，为结构设计提供了依据。

目前风吹雪数值模拟多为对简单屋面积雪分布的分析，缺少其在大型场馆上的应用。本文基于Euler-Euler方法，以某大型体育场为研究对象，研究不同风速、风向下屋面各分区雪荷载的分布规律。数值模拟时采用控制方程计算风场信息及雪相分布，通过软件FLUENT中的UDF接口获取壁面处的摩擦速度、雪体积分数等变量，并根据积雪预测模型判断壁面处雪颗粒的沉积、侵蚀，最后经过源项反馈到控制方程中。为验证风吹雪模拟方法的有效性和实用性，对球壳屋面雪荷载分布进行模拟，与已有文献结果进行对比，分析了积雪分布形式和积雪分布系数的分布情况。

1 风吹雪的数值模拟方法

1.1 控制方程

本文采用计算流体力学软件FLUENT中的Mixture多相流模型进行模拟计算，选取Realizable k-\epsilon湍流模型，双精度求解。

在Mixture多相流模型中空气相和雪相共用一套连续方程和动量方程，方程中通过相对漂移速度考虑了雪相对空气流动的影响，并通过弥散相连续方程计算流场中雪颗粒的体积分数。

混合相连续方程如下所示：

\frac{\partial }{\partial t}({\rho }_{\mathrm{m}})+\nabla ·({\rho }_{\mathrm{m}}{\mathit{v}}_{\mathrm{m}})=\mathrm{0}（1）

{\rho }_{\mathrm{m}}=f_{\mathrm{a}}{\rho }_{\mathrm{a}}+f_{\mathrm{s}}{\rho }_{\mathrm{s}}（2）

{\mathit{v}}_{\mathrm{m}}=\frac{\mathrm{1}}{{\rho }_{\mathrm{m}}}\left(f_{\mathrm{a}}{\rho }_{\mathrm{a}}{\mathit{v}}_{\mathrm{a}}+f_{\mathrm{s}}{\rho }_{\mathrm{s}}{\mathit{v}}_{\mathrm{s}}\right)（3）

式中：f为体积分数；\rho为密度；\mathit{v}为速度向量；\nabla为哈密顿算子，\nabla =\mathit{i}\frac{\partial }{\partial x}+\mathit{j}\frac{\partial }{\partial y}+\mathit{k}\frac{\partial }{\partial z}，\mathit{i}、\mathit{j}、\mathit{k}为单位向量；下标a、s、m分别表示空气相、雪相、混合相。

混合相动量方程如下所示：

\begin{array}{l}\frac{\partial ({\rho }_{\mathrm{m}}{\mathit{v}}_{\mathrm{m}})}{\partial t}+\nabla ·({\rho }_{\mathrm{m}}{\mathit{v}}_{\mathrm{m}}{\mathit{v}}_{\mathrm{m}})=-\nabla p+\nabla ·\left[{\mu }_{\mathrm{m}}(\nabla {\mathit{v}}_{\mathrm{m}}+\nabla {\mathit{v}}_{\mathrm{m}}^{\mathrm{T}})\right]\\ +{\rho }_{\mathrm{m}}\mathit{g}+\mathit{F}+\nabla ·\left({\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{2}}}{f}_k{\rho }_k\mathit{v}{}_{dr,k}{}^{}\mathit{v}{}_{dr,k}{}^{}\right)\end{array}（4）

式中：p为流场压力；{\mu }_{\mathrm{m}}为空气相和雪相的混合动力粘度；g为重力加速度；v_{\mathrm{d}\mathrm{r},\mathrm{k}}=v_{\mathrm{k}}-v_{\mathrm{m}}，为相对漂移速度。整体形式上，相较单相动量方程，混合相动量方程额外引入扩散项（方程右侧最后一项）来考虑相间滑移导致的动量扩散。

雪相连续方程如下所示：

\frac{\partial (f_{\mathrm{s}}{\rho }_{\mathrm{s}})}{\partial t}+\nabla ·(f_{\mathrm{s}}{\rho }_{\mathrm{s}}{\mathit{v}}_{\mathrm{m}})=-\nabla ·(f_{\mathrm{s}}{\rho }_{\mathrm{s}}{\mathit{v}}_{\mathrm{d}\mathrm{r},\mathrm{s}})+S_{\mathrm{q}}（5）

式中：S_{\mathrm{q}}为与雪颗粒侵蚀、沉积相关的源项，壁面处的积雪通量变化通过S_{\mathrm{q}}考虑，具体计算方法见1.2节。本文不考虑雪相的升华等与空气相成分互相转换的情况。

1.2 积雪预测模型

积雪预测模型是风致雪漂移正确模拟的重要因素之一，其决定了壁面处雪颗粒的沉积与侵蚀。本文采用NAAIM等^［12］提出的基于壁面剪切关系建立的积雪预测模型，其考虑了摩擦速度和雪质量浓度的双重影响，避免了对壁面积雪通量的过高估计。单位时间内壁面上雪颗粒侵蚀（沉积）质量通量q_{\mathrm{e}\mathrm{r}\mathrm{o}}（q_{\mathrm{d}\mathrm{e}\mathrm{p}}）的计算公式如下：

q_{\mathrm{e}\mathrm{r}\mathrm{o}}=A_{\mathrm{e}\mathrm{r}\mathrm{o}}(u_{\mathrm{*}}^{\mathrm{2}}-u_{\mathrm{*}\mathrm{t}}^{\mathrm{2}}),u_{\mathrm{*}}>u_{\mathrm{*}\mathrm{t}}（6）

（7）

式中：A_{\mathrm{e}\mathrm{r}\mathrm{o}}为常数，与颗粒间黏合程度有关，取7×10^-4 kg·m^-4·s；u_{\mathrm{*}}=\sqrt[]{\tau /{\rho }_{\mathrm{a}}}，为壁面处摩擦速度，式中τ为壁面剪切应力；u_{\mathrm{*}\mathrm{t}}为阈值摩擦速度；C为雪质量浓度；w_{\mathrm{f}}为雪的沉降速度。模拟计算中u_{\mathrm{*}\mathrm{t}}取0.25 m·s^-1，{\rho }_{\mathrm{s}}为雪颗粒密度取250 \mathrm{k}\mathrm{g}\cdot {\mathrm{m}}^{-\mathrm{3}}，w_{\mathrm{f}}取0.2 m·s^-1。

根据公式可知，壁面上雪颗粒的侵蚀或沉积是由壁面上的摩擦速度决定的，当u_{\mathrm{*}}>u_{\mathrm{*}\mathrm{t}}时，壁面上积雪发生侵蚀，进入计算域；当时，壁面上积雪发生沉积，离开计算域；当u_{\mathrm{*}}=u_{\mathrm{*}\mathrm{t}}时积雪面保持平衡，壁面雪通量为零。

近壁面网格中雪相连续方程的源项S_{\mathrm{q}}的表达式为：

S_{\mathrm{q}}=\frac{q_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}}{h_{\mathrm{p}}}（8）

q_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}=q_{\mathrm{e}\mathrm{r}\mathrm{o}}+q_{\mathrm{d}\mathrm{e}\mathrm{p}}（9）

式中：h_{\mathrm{p}}为壁面第一层网格高度。

1.3 边界条件

计算时需给定合理的边界条件，才能得出所需流场的解。计算域各边界条件如表1所示，其中计算域入口选用速度入口（velocity-inlet），入流风速采用指数型风剖面计算，数值模拟过程中分别考虑西安日平均风速以及极端情况西安50年重现期风速（最大降雪量发生的同时出现可能的最高风速），即参考高度z_{\mathrm{0}}=10 m处的风速u_{\mathrm{0}}分别取3.0 m·s^-1和23.6 m·s^-1，入口湍流采用文献［13］中的经验公式计算。

考虑到钝体绕流的分离与附着现象，局部风雪流难以达到饱和，吹雪分层界面及各层质量交换情况比较复杂^［14］，故不在入流边界处对跃移或悬移进行区分，入流雪体积分数取常数。

计算域出口采用pressure-outlet，以防止出口回流对结果的影响；计算域顶面及侧面采用对称边界条件symmetry；底面采用无滑移壁面wall，并使用壁面函数计算近壁面湍流。

2 计算模型及数据处理

2.1 几何建模及网络划分

该大型场馆的屋面由网壳屋面（大屋面）以及索膜屋面（小屋面）两部分组成（图2），大小屋面均为双曲型；屋面波谷处沿x、y方向屋面斜率比较平缓，波峰处较陡峭；小屋面为单向波纹状，波纹与y向平行，沿x向形成波形起伏；大屋面内边缘与小屋面外边缘有4 m左右的高差。

计算模型根据实际尺寸进行1∶1的建模。计算域的长、宽、高分别为4 500 m×2 000 m×500 m。考虑计算速度及结果的准确性，采用区域分块技术进行网格划分，将计算区域分为建筑周边区域和其他区域，建筑周边区域采用非结构化网格，其他区域采用结构化网格，并对建筑周边区域网格进行加密。计算域最小网格尺寸为0.2 m，网格总量为580万，模型阻塞比小于5%，建筑表面网格划分情况如图3所示。

考虑建筑的对称性，仅模拟风向角的范围在0°~180°之间，每30°一个工况（图4），根据西安风玫瑰图，西安主导风向为东北方向，第二主导风向为西南方向，故额外选定45°、135°方向进行模拟计算，考虑0°、30°、45°、60°、90°、120°、135°、150°、180°共9个风向角。

将屋面划分为232个分区，大屋面中160个，小屋面中72个，分区如图4所示。选取屋面典型位置处的12个分区，编号分别为X-1~X-12，图中蓝色椭圆内部为小屋面，外部为大屋面。

2.2 数据处理

荷载规范中提出建筑结构屋面上的积雪荷载由当地的基本雪压s_{\mathrm{0}}和屋面上的积雪分布系数{\mu }_{\mathrm{r}}决定，{\mu }_{\mathrm{r}}等于屋面水平投影面积上的雪荷载s_{\mathrm{h}}与基本雪压s_{\mathrm{0}}的比值，即基本雪压转换成屋面积雪荷载的换算系数。故考虑积雪侵蚀沉积的屋面积雪分布系数由下式计算：

{\mu }_{\mathrm{r}}=\frac{s_{\mathrm{h}}}{s_{\mathrm{0}}}=\frac{\mathrm{1}}{\mathrm{c}\mathrm{o}\mathrm{s}\alpha }\left(\mathrm{c}\mathrm{o}\mathrm{s}\alpha +\frac{\mathrm{\Delta }h{\rho }_{\mathrm{s}}g}{\mathrm{1}\mathrm{000}{s}_{\mathrm{0}}}\right)（10）

式中：\alpha为屋面切线与水平面的夹角；\mathrm{\Delta }h=q_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}\mathrm{\Delta }t/{\rho }_{\mathrm{s}}，为\mathrm{\Delta }t时间内屋面积雪迁移引起的雪高度变化值。由公式可以看出{\mu }_{\mathrm{r}}大于1时表示雪颗粒沉积，小于1时表示积雪侵蚀，与荷载规范中其代表的意义一致。

定义屋面沉积率^［9］为{\eta }_{\mathrm{d}\mathrm{e}\mathrm{p}}，用以定量表述屋面表面雪颗粒沉积情况，公式如下所示：

{\eta }_{\mathrm{d}\mathrm{e}\mathrm{p}}=\frac{A_{\mathrm{d}\mathrm{e}\mathrm{p}}}{A_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}}\times \mathrm{100}\mathrm{\%}（11）

式中：A_{\mathrm{d}\mathrm{e}\mathrm{p}}为屋面雪颗粒沉积区域的面积；A_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}为屋面表面的总面积。

此外为评估屋面雪荷载的不均匀性，参考文献［15］中评估双坡屋面雪荷载分布不均匀性的方法，定义大跨屋面雪荷载不均匀性系数C_{\mathrm{u}\mathrm{n}}为：

C_{\mathrm{u}\mathrm{n}}=\frac{\left|{\overline{s}}_{{}_{\mathrm{e}\mathrm{r}\mathrm{o}}}-{\overline{s}}_{{}_{\mathrm{d}\mathrm{e}\mathrm{p}}}\right|}{s_{\mathrm{0}}}（12）

式中：{\overline{s}}_{{}_{\mathrm{e}\mathrm{r}\mathrm{o}}}为侵蚀区域的平均雪压；{\overline{s}}_{{}_{\mathrm{d}\mathrm{e}\mathrm{p}}}为沉积区域的平均雪压。

3 数值模拟方法合理性验证

为验证数值模拟方法的合理性，对与文献［16］中球壳模型大小一致的几何模型进行数值模拟分析，并与文献中的风洞试验结果进行对比。

该球壳模型屋檐高H=0.2 m，跨度为5H，矢跨比为1/5。计算域尺寸为155H×50H×10H，球壳屋面位于计算域1/3处。模拟中地貌为B类，参考高度为1 m处的风速，为3.0 m·s^-1，雪密度取250 \mathrm{k}\mathrm{g}\cdot {\mathrm{m}}^{-\mathrm{3}}，雪颗粒直径取0.000 15 m，阈值摩擦速度取0.3 m·s^-1。

图5为模拟与试验无量纲雪深分布对比图，图5a）为本文模拟结果，图5b）为文献［16］风洞试验结果。由图可知，数值模拟得出的球壳屋面积雪分布形式与风洞试验结果一致：积雪在前部、后部堆积，中部处于被侵蚀区域没有积雪，由于颗粒的休止角影响，屋面边缘处积雪厚度均不大。提取模型中轴线上的积雪分布绘于图6，图中横坐标为屋面上测点位置与屋檐高度的比值x/H；纵坐标为无量纲雪深S/H（实际雪深与屋檐高度的比值），黑色曲线为屋面曲线。相比文献［16］，屋面迎风侧积雪的模拟结果位置更靠前，且积雪深度小于文献中结果；屋面背风侧积雪的覆盖面积大于文献结果，两者最大雪深相差7%。

综上所述，数值模拟的球壳整体积雪分布形式及积雪厚度均与风洞试验结果基本一致，验证了本文的数值模拟方法在模拟曲面的积雪分布时是可行的。

4 计算结果及分析

为方便描述，文中日平均风速下某风向角工况用AV-xxx表示，50年重现期风速下某风向角工况用BA-xxx表示。

4.1 典型风向下摩擦速度分析

限于篇幅本文仅给出风向角为0°、45°、90°时，日平均风速及50年重现期风速作用下的屋面u_{\mathrm{*}}的分布情况（图7）。

由图7可以看出，屋面u_{\mathrm{*}}的分布规律与来流风向密切相关，风速对其影响有限。风速不变时，来流风两侧的屋脊处的u_{\mathrm{*}}远大于u_{\mathrm{*}\mathrm{t}}，该区域积雪将被侵蚀，如当风向角为45°时大屋面的屋脊B、D处；迎风凸起处的u_{\mathrm{*}}处于中等水平，与u_{\mathrm{*}\mathrm{t}}相差不大，该区域雪颗粒是否被侵蚀，与来流风速大小相关，如当风向角为45°时的屋脊A和风向角为90°时的波峰A处；屋面背风侧u_{\mathrm{*}}普遍小于u_{\mathrm{*}\mathrm{t}}，易出现积雪堆积现象，如当风向角为0°时的波谷B、风向角为45°时的屋脊C以及风向角为90°时的波峰B；波谷处于上游时其u_{\mathrm{*}}仍处于较低值，易出现积雪堆积现象，如当风向角为0°时的波谷A。

来流风经波谷穿过时（风向角为0°），大屋面对小屋面的遮挡较弱，此时小屋面内边缘迎风位置（即气流分离处）u_{\mathrm{*}}较大，紧邻大屋面的外边缘u_{\mathrm{*}}较小，小屋面下游u_{\mathrm{*}}更小；当风向角为90°时，来流风与小屋面波纹垂直且受波峰A阻挡，小屋面上游u_{\mathrm{*}}比下游u_{\mathrm{*}}小；当风向角为45°时，小屋面u_{\mathrm{*}}的分布规律接近于风向角为0°时小屋面u_{\mathrm{*}}的分布规律。

50年重现期风速下的屋面u_{\mathrm{*}}的分布形式与日平均风速时基本相同，屋面u_{\mathrm{*}}值的变化受风向影响。随风速增大，屋面迎风区的u_{\mathrm{*}}逐渐增大超过u_{\mathrm{*}\mathrm{t}}，雪颗粒开始被侵蚀；屋面背风区的u_{\mathrm{*}}仍然较小，雪颗粒持续在背风侧堆积。

4.2 屋面积雪分布系数

选取西安地区基本雪压0.25 kN·m^-2作为基准值，计算风吹雪8 h后屋面的积雪分布系数。积雪分布系数为无量纲值，计算结果如图8所示，图中积雪分布系数大于等于1.00的显示为红色。

由图可知，日平均风速下，屋面积雪分布系数多处于1.00±0.10范围内，积雪分布相对均匀，不同风向下积雪的沉积、侵蚀与屋面摩擦速度分布密切相关。

工况AV-000时，来流风正对波谷A，与小屋面波纹平行，波谷A、B位置处雪颗粒发生沉积，积雪分布系数最大值为1.13，位于下游波谷B；受空气绕流影响，屋面两侧风速较大积雪受到侵蚀，积雪分布系数多在0.95附近；大屋面与小屋面相接处高差4 m，形成高低屋面，使得两屋面相接处附近风速较低，积雪沉积严重，积雪分布系数最大为1.32。工况AV-090时，来流风正对波峰A，来流两侧的波谷A、B处积雪被侵蚀，由于下部结构南北不对称，所以波谷A、B的侵蚀状况不同。与工况AV-000相比，工况AV-090下的迎风面更陡，屋檐和屋顶处气流分离更加剧烈，使得大屋面波峰A的屋檐和屋顶处积雪出现侵蚀现象，最大积雪分布系数出现在下游波峰近屋脊位置，值为1.48；此时风向与小屋面波纹垂直，波纹对雪颗粒滑动的阻挡达到最大化，小屋面雪颗粒沉积区域增大，最大积雪分布系数为1.03。工况为AV-045时，来流风接近屋脊A方向，屋脊C位于背风侧，大屋面屋脊B、D处积雪被侵蚀，最大侵蚀量为17%，雪颗粒堆积集中出现在屋脊C处，该位置最大积雪分布系数为1.08；此时风向与小屋面波纹斜交，小屋面雪荷载沉积区域的范围介于在工况AV-000与工况AV-090下小屋面雪荷载沉积区域之间。

随风速增大，积雪被侵蚀的面积逐渐扩大，积雪分布系数变化梯度增大，积雪不均匀分布加剧。50年重现期风速下，工况为BA-000时，大屋面上游波谷处的积雪被侵蚀，雪颗粒的沉积集中在背风侧波谷处，大屋面最大的积雪分布系数为2.99，而屋面两侧波峰的积雪侵蚀量达80%；相比大屋面，小屋面波峰处积雪侵蚀量变化较小，但下游波谷处积雪明显增加。工况为BA-045时，大屋面仅个别区域出现雪颗粒沉积现象，最大雪压出现在波谷B处，积雪分布系数为1.61；此时小屋面波纹无法阻挡风对积雪的侵蚀，雪颗粒仅在波谷B附近产生沉积，最大积雪分布系数为2.15。工况为BA-090时，大屋面雪颗粒沉积集中在下游波峰处，积雪分布系数在1.10附近；迎风两侧波谷侵蚀量接近，均在50%左右，小屋面在上下游波峰处均出现雪颗粒堆积现象，积雪分布系数值最大为1.36，出现在下游内圈附近。

图9为风向角为0°、90°时，体育场迎风断面处无量纲风速分布以及风速矢量。由图可以看出，风向角为0°时，来流风与小屋面波纹平行，气流从波纹中流过未出现回流，由于大屋面的遮挡，小屋面与大屋面相接处风速几乎为0，大屋面上方积雪被吹落时，将在相接处产生堆积。风向角为90°时，小屋面高度低于大屋面高度，其上产生明显的回流旋涡，且由于来流风与小屋面波纹垂直，小屋面上方风场较复杂，波纹底部风速接近于0，此处产生积雪堆积。

不同风向角下，屋面沉积率及不均匀系数如图10所示。日平均风速下，屋面沉积率最大值为47.8%，出现在工况AV-090下，屋面不均匀性系数最大为0.16，出现在工况AV-135下。随着风速增加，屋面沉积率减小，不均匀系数增大。50年重现期风速下，最大沉积率值为12.8%，出现在工况BA-120下；屋面不均匀性系数最大值为1.16，出现在工况BA-000下。与日平均风速相比，50年重现期风速下，屋面积雪沉积率小于日平均风速时屋面积雪沉积率的1/2，但局部积雪沉积量更大，这是因为上游更多的雪被吹向摩擦速度较低的区域，由此造成在高风速下屋面积雪不均匀分布更加明显，50年重现期风速下屋面积雪不均匀系数比日平均风速下屋面积雪不均匀系数大一个数量级。

考虑到不同风速下屋面积雪分布情况存在显著差异，尽管50年重现期基本雪压和50年重现期基本风速同时出现的概率很小，为确保结构安全，结构设计时也应综合考虑常遇风和基本风速下屋面雪分布的影响。

4.3 积雪分布系数随风向的变化

对屋面积雪分布系数进行统计，计算屋面积雪分布系数的平均值，屋面积雪分布系数最大值与平均值随风向角的变化如图11所示。由图可知：不同风速下，屋面积雪分布系数随风向角变化不同；日平均风速下，积雪分布系数最大值为1.48，出现在风向角为90°时，积雪分布系数平均值随风向角的变化不大，处于0.96~1.00范围内；50年重现期风速下，积雪分布系数最大值随风向角增大呈先减小后增大趋势，风向角为90°时为转折点，当风向角为180°时积雪分布系数最大，其值为4.70；风速增大时屋面积雪侵蚀量增大，50年重现期风速下积雪分布系数平均值在0.40~0.70之间，为日平均风速时积雪分布系数平均值的2/3左右。

在屋面上选取12个典型位置（图4），进一步分析不同位置处积雪分布系数与风向角的关系。其中分区X-1~X-8位于大屋面，分区X-9~X-12位于小屋面；分区X-1、X-5、X-9、X-11位于屋面波谷处，分区X-3、X-7、X-10、X-12位于屋面波峰处；分区X-2、X-4、X-6、X-8位于屋脊。不同风速下典型位置处积雪分布系数随风向变化曲线如图12、图13所示。

对比图12a）、图13a）可以看出：较低风速下，分区X-1和X-5分别在风向角为30°~120°和风向角为90°~135°时，积雪处于被侵蚀状态；位于小屋面的分区X-9和X-11在多数风向角下出现雪颗粒沉积现象，但沉积量不大。在高风速下，积雪分布系数规律性更明显，大小屋面上位于相同波谷区域的积雪分布系数随风向的变化趋势相似；波谷受到侵蚀的风向角范围更大，分区X-1和X-5分别在风向角为0°~120°时、风向角为60°~180°时处于被侵蚀状态，位于小屋面的分区X-9和X-11仅在个别风向下出现雪颗粒沉积现象；波谷处于背风侧时，可能会出现较大的沉积量，分区X-1在风向角为135°时积雪分布系数为2.49，是日平均风速下积雪分布系数的2.37倍，分区X-11在风向角为45°时积雪分布系数为2.15，是日平均风速下积雪分布系数的2.13倍。

由图12b）和图13b）可知，波峰位于背风侧时易出现雪颗粒沉积现象。日平均风速和50年重现期风速下的分区X-3和X-7均在风向角为90°时出现雪颗粒沉积现象，风速增大后，波峰处积雪受侵蚀更加显著，当波峰位于侧风面时，积雪分布系数在0附近；由于大屋面对小屋面的遮挡，低风速下小屋面上的分区X-10、X-12仅在个别风向角时被侵蚀，在高风速时该现象减弱，分区X-10在风向角为60°~120°时雪颗粒出现沉积现象，而分区X-12始终被侵蚀。

由图12c）和图13c）可知，屋脊位置积雪分布系数变化幅度较大。屋脊A、B、C、D随风向变化轮番出现雪颗粒侵蚀与沉积现象，高风速下此现象更加明显，这与四个屋脊的位置有关。比如在风向角为60°时，分区X-4与X-8分别位于波谷的上下游，此时分区X-4处雪颗粒堆积而分区X-8处积雪被侵蚀；分区X-6在风向角为30°~90°时以及分区X-8在风向角为120°~150°时出现雪颗粒堆积现象。

由以上结果可以看出，高风速下，积雪分布系数对风向更加敏感，50年重现期风速下，大屋面多数区域处的积雪分布系数随风向的变化趋势与日平均风速下积雪分布系数的变化趋势类似，积雪分布系数随风向角变化起伏较大，屋面雪颗粒沉积量和积雪侵蚀量均较大。随着风速的增加，小屋面的波纹形状对雪颗粒的侵蚀沉积影响减弱，其上游区域积雪被侵蚀，其积雪分布系数随风向角的变化与大屋面积雪分布系数变化一致。

5 结论

本文基于Mixture多相流模型，研究了不同风速、风向下某大型屋面积雪分布情况，得到以下结论：

（1） 大跨屋面的雪荷载分布受风速、风向影响较大，荷载规范中没有考虑相关影响，结构设计时有必要综合考虑不同风速、风向对雪荷载的影响。

（2） 相同风速下，屋面背风侧雪颗粒发生沉积，侧风面积雪受到侵蚀，在低处屋面容易产生积雪堆积。低风速时，迎风侧可能会出现雪颗粒堆积现象，屋面局部造型对积雪分布影响较大。高风速时，风向对积雪分布影响更明显，局部沉积量增大。

（3） 随着风速增加，屋面沉积率减小，不均匀系数增大。在日平均风速下，从平均意义上来说，屋面雪颗粒发生沉积的面积是50年重现期风速下屋面雪颗粒发生沉积面积的3.6倍，此时屋面的积雪不均匀系数比日平均风速下的积雪不均匀系数大一个数量级。

（4） 在日平均风速下，积雪分布系数平均值在0.96~1.00之间；在50年重现期风速下积雪分布系数平均值在0.40~0.70之间，相比日平均风速作用下，50年重现期风速下积雪总量减少1/3，但最不利雪压增大两倍以上。

（5） 目前大跨屋盖风吹雪模拟中未考虑雪融化以及二次堆积造成的影响，后期研究中可以增加对相变过程的模拟，实现对积雪的长期预测，为结构全生命周期防灾提供依据。

参考文献

[1] 徐晓明， 陈伟， 史炜洲， 西安国际足球中心屋盖钢结构设计［J］. 建筑结构， 2023， 53（1）： 47-53. DOI：10.19701/j.jzjg.SHY2326.

[2] XU Xiaoming， CHEN Wei， SHI Weizhou， et al. Structural design on steel roof of Xi’an international football center［J］. Building Structure， 2023， 53（1）： 47-53. DOI：10.19701/j.jzjg.SHY2326.（in Chinese）

[3] GB 50009—2012 建筑结构荷载规范［S］．

[4] GB 50009—2012 Load code for the design of building structures［S］. （in Chinese）

[5] BEYERS J H M， SUNDSB P A， HARMS T M. Numerical simulation of three-dimensional， transient snow drifting around a cube［J］. Journal of Wind Engineering and Industrial Aerodynamics， 2004， 92（9）： 725-747. DOI：10.1016/j.jweia.2004.03.011.

[6] 孙晓颖， 洪财滨， 武岳， 建筑物周边风致雪漂移的数值模拟研究［J］. 工程力学， 2014， 31（4）： 141-146.

[7] SUN Xiaoying， HONG Caibin， WU Yue， et al. Numerical simulation of snow drifting around building model［J］. Engineering Mechanics， 2014， 31（4）： 141-146. （in Chinese）

[8] 孙晓颖， 洪财滨， 武岳. 典型形式大跨度屋盖风雪漂移的数值模拟［J］. 振动与冲击， 2014， 33（18）： 36-42. DOI：10.13465/j.cnki.jvs.2014.18.007.

[9] SUN Xiaoying， HONG Caibin， WU Yue. Numerical simulation of snow drifting on typical long-span roofs［J］. Journal of Vibration and Shock， 2014， 33（18）： 36-42. DOI：10.13465/j.cnki.jvs.2014.18.007.（in Chinese）

[10] 刘博雅， 张清文， 刘盟盟， 高跨为双坡屋面的高低屋面低跨平屋面积雪分布研究［J］. 建筑结构学报， 2019， 40（6）： 48-55. DOI：10.14006/j.jzjgxb.2018.0370.

[11] LIU Boya， ZHANG Qingwen， LIU Mengmeng， et al. Research of snow distribution on lower roof of step roofs with gable-roof on high part［J］. Journal of Building Structures， 2019， 40（6）： 48-55. DOI：10.14006/j.jzjgxb.2018.0370.（in Chinese）

[12] ZHANG Z B， MA W Y， LI Q， et al. Snow load shape coefficients and snow prevention method for stepped flat roofs［J］. Applied Sciences， 2023， 13（22）： 12109. DOI：10.3390/app132212109.

[13] 孙芳锦， 孙红玉， 于珊珊， 风雪共同作用下大跨度双曲屋盖数值模拟研究［J］. 应用力学学报， 2023， 40（4）： 926-933. DOI：10.11776/j.issn.1000-4939.2023.04.023.

[14] SUN Fangjin， SUN Hongyu， YU Shanshan， et al. Numerical simulation of large-span hyperbolic roofs under the combined action of wind and snow［J］. Chinese Journal of Applied Mechanics， 2023， 40（4）： 926-933. DOI：10.11776/j.issn.1000-4939.2023.04.023.（in Chinese）

[15] 周晅毅， 顾明， 李雪峰. 大跨度屋盖表面风致雪压分布规律研究［J］. 建筑结构学报， 2008， 29（2）： 7-12. DOI：10.14006/j.jzjgxb.2008.02.002.

[16] ZHOU Xuanyi， GU Ming， LI Xuefeng. Study on wind-induced snow pressure on large-span roof［J］. Journal of Building Structures， 2008， 29（2）： 7-12. DOI：10.14006/j.jzjgxb.2008.02.002.（in Chinese）

[17] 徐爱民， 陈志强， 马永兴， 长沙机场T3航站楼雪荷载数值模拟及取值［J］. 重庆大学学报， 2022， 45（12）： 116-124. DOI：10.11835/j.issn.1000-582X.2021.059.

[18] XU Aimin， CHEN Zhiqiang， MA Yongxing， et al. Numerical simulation of snow loading on Changsha airport terminal 3［J］. Journal of Chongqing University， 2022， 45（12）： 116-124. DOI：10.11835/j.issn.1000-582X.2021.059.（in Chinese）

[19] 卢春玲， 陈锦焜， 刘宇杰. 大跨度屋盖风致积雪重分布的数值模拟［J］. 建筑结构， 2021， 51（15）： 107-114. DOI：10.19701/j.jzjg.2021.15.017.

[20] LU Chunling， CHEN Jinkun， LIU Yujie. Numerical simulation of wind-induced snow redistribution on large-span roof［J］. Building Structure， 2021， 51（15）： 107-114. DOI：10.19701/j.jzjg.2021.15.017.（in Chinese）

[21] NAAIM M， NAAIM-BOUVET F， MARTINEZ H. Numerical simulation of drifting snow： Erosion and deposition models［J］. Annals of Glaciology， 1998， 26： 191-196. DOI：10.3189/1998aog26-1-191-196.

[22] 王福军. 计算流体动力学分析： CFD软件原理与应用［M］. 北京： 清华大学出版社， 2004.

[23] WANG Fujun. Computational fluid dynamics analysis： Principle and application of CFD software［M］. Beijing： Tsinghua University Press， 2004. （in Chinese）

[24] 刘多特， 李永乐， 汪斌. 风雪绕流数值模拟的积雪预测模型研究［J］. 工程力学， 2016， 33（8）： 122-131. DOI：10.6052/j.issn.1000-4750.2015.01.0019.

[25] LIU Duote， LI Yongle， WANG Bin. A numerical prediction model for snow accumulation caused by ambient snowdrift［J］. Engineering Mechanics， 2016， 33（8）： 122-131. DOI：10.6052/j.issn.1000-4750.2015.01.0019. （in Chinese）

[26] 刘志祥， 余志祥， 何欢， 典型单跨低坡度双坡屋面的风致雪堆试验研究［J］. 建筑结构学报， 2024， 45（5）： 136-144. DOI：10.14006/j.jzjgxb.2023.0094.

[27] LIU Zhixiang， YU Zhixiang， HE Huan， et al. Experimental studies of wind-induced snowdrifts on typical single-span low-slope gable roofs［J］. Journal of Building Structures， 2024， 45（5）： 136-144. DOI：10.14006/j.jzjgxb.2023.0094.（in Chinese）

[28] 王建烁， 刘红波， 陈志华， 球形屋顶积雪分布风洞试验研究［J］. 天津大学学报（自然科学与工程技术版）， 2021， 54（7）： 711-718. DOI：10.11784/tdxbz202001020.

[29] WANG Jianshuo， LIU Hongbo， CHEN Zhihua， et al. Wind tunnel test of snow distribution on spherical roofs［J］. Journal of Tianjin University （Science and Technology）， 2021， 54（7）： 711-718. DOI：10.11784/tdxbz202001020.（in Chinese）