多孔介质方法对铁路沿线沙障模拟的适用性分析

王 连，程建军，智凌岩，辛国伟

(石河子大学水利建筑工程学院，新疆石河子 832003)



多孔介质方法对铁路沿线沙障模拟的适用性分析

王 连，程建军，智凌岩，辛国伟

(石河子大学水利建筑工程学院，新疆石河子 832003)

为模拟细孔透隙式铁路沙障的流场特征，鉴于此类模拟建模的困难性，采用CFD数值模拟程序，基于三维模拟的多孔介质条件对阻沙固沙网和防风沙直立栅栏的流场特征及压力变化进行模拟分析。通过边界条件和参数的合理控制得到：直立栅栏流场的数值模拟结果与风洞实验下的流场速度变化相似；由于栅栏的整流作用，在多孔介质模拟方法下沙障后出现明显的减速恢复区，并在其后逐渐恢复为与入口相同的风速廓线流；结合全断面PE网风洞试验，采用多孔介质模型可以得到与试验数据相吻合的模拟结果，且压降随入口速度的增加明显增大；由模拟结果可知，通过参数的合理控制，多孔介质方法可较好地反映此类沙障的流场变化。

多孔介质；沙障；压降；流场

1 概述

在我国尤其是西北地区，荒漠化严重，风沙灾害频繁，早在20世纪60年代，我国就开始风沙灾害的治理工作，从早期的植物防护到如今的机械防护，风沙治理工作取得了相对理想的进展。但面对日渐恶劣的环境加之不同季节特殊地理位置的气候条件，科学工作者仍需针对防风沙障碍物作用效果做进一步研究。目前，PE网应用较为广泛，其防风固沙性能优异，且具有优异的抗老化性能，连片使用可显著增大地面粗糙度，降低来流风速；独特的上疏下密结构可引导风向沙障较梳的中上部通过，显著减弱风对地表的挖蚀作用，具有减风、阻沙、导风的综合效果[1]。尼龙阻沙网作为一种新型的防沙材料，不仅具有防风阻沙效果，还具有疏导沙的作用，在防沙治沙工程中具有广泛应用前景。

阻沙栅栏又称高立式沙障，在我国风沙地区广泛使用，直立栅栏通过改变局部流场，使一定范围内的近地表风速明显降低，使风沙流携沙于此堆积[2]。栅栏防护体系的最主要作用是降低风速，增加风沙流饱和度，加速栅栏周围沙粒的沉积，抑制风沙灾害发生，其防护效益通常用削减风速的大小来评价[3]。防风沙栅栏主要用于铁路沿线对风沙灾害的防治，设置在防护体系的最外侧，与铁路走向平行，以减少风沙运动对铁路沿线运输的影响，PE网和防风沙栅栏均是沙漠地区阻沙固沙的理想新型材料。故本文通过计算流体动力学(CFD)数值模拟软件对二者进行数值模拟，减少实际测试或风洞试验获取数据过程中的人力、物力投入。经过多次模拟对比分析，针对不同模型寻找最佳模拟方法，通过对模拟计算过程中不同阶段参数的合理控制，并结合风洞实验的数据加以修正，最终得到合理的计算方法和湍流模型，并分析障后流场的变化和通过防风沙障碍物前后压力的变化，凸显防风沙障碍物以降低风速和减少湍流强度为目标的防治原理，为从事风沙研究的工作者提供一定的参考。

2 试验设计分析

2.1 风洞试验设计

试验在中国科学院新疆生态与地理研究所风洞实验室进行，该风洞全长16.2 m，由动力段、试验段和扩散段构成，试验段长8 m，横截面为宽1.3 m、高1 m的矩形，属直流吹气式风洞。将测试装置固定于风洞底板，皮托管伸入孔洞，调整至合适位置并固定，将风洞入口风速调整至不同待测值，待皮托管读数仪上显示的压强数值稳定后记录数值，依次测试各预留孔洞处的压强数值，即可得到相应风速下PE网前后各处的压差值。试验用PE网为甘肃金海阻沙固沙新材料有限公司生产的JGSPE200He18G型固沙网，其孔隙率约60%，孔规格0.79 mm×0.83 mm。风洞试验布置如图1所示。

图1 风洞试验布置示意

2.2 全断面PE网数值模拟

由于实体模型对划分网格带来诸多不便，故对于三维的PE网，采用数值模拟软件自带的多孔介质模型模拟其压降随入口速度的变化，并与试验数据进行对比分析。

2.2.1 多孔介质模型

多孔介质通常具有30%～60%孔隙度，模拟方式具有多样性。CFD中多孔介质模型是在定义多孔介质的区域结合了一个根据经验假设为主的流动阻力，即在动量方程上增加了一个动量源项，其中源项由两部分组成，黏性损失项和内部损失项。1856年提出的Darcy关系式第一次描述了多孔介质内单向流速与压降之间的关系式[4]

(1)

式中，Δp为单相流压降，kPa；L为多孔介质区域长度，m；μ为运动黏度，Pa·s；k为渗透率，m2；u为单相流表观速度，m·s-1。

Forchheimer指出，多孔介质中单项流的压降实际上由黏性力和惯性力两部分组成，低雷诺数下惯性力的作用可忽略不计，而高雷诺数下则影响显著；Ergun总结前人的研究成果，认为高雷诺数下多孔介质中单向流的压降应表述为速度一次项和二次项之和的形式，并提出Ergun关系式[5-8]。

2.2.2 模型的建立

数值模拟几何模型如图2所示。

图2 数值模拟几何模型

模型长度方向取41.0 mm，宽度方向取11.4 mm，高度方向取9.8 mm，多孔区域位于长度中线位置，Δn=0.27 mm，采用全断面布设多孔区域以验证试验结果。

2.2.3 网格划分

利用计算流体动力学模拟软件对既有几何模型的网格进行划分，其中对于多孔介质区域，将其独立于计算区域进行网格划分，并对靠近多孔介质壁面区域进行网格加密处理，统一采用六面体网格即Hexa Meshing方式进行划分，最终得到286 110个节点、830 900个网格面和272 500个网格单元，网格最小正交质量为1.0。划分网格结果如图3所示。

图3 模型网格划分

2.2.4 模拟计算及边界条件

多孔介质条件亦称分布阻力法，是将流动区域中固体结构的作用看作是附加在流体上的分布阻力，此时多孔模型的阻碍作用不复存在，只需设置相关参数，而数值模拟软件默认多孔介质区域的速度为依据体积流量来推断的表面速度，同时也可以赋予多孔介质内部速度值。多孔介质对湍流的影响是近似的流体力学湍流动量方程

(2)

方程中第一项为时间项；第二项为对流项；第三项为扩散项；最后一项为源项，方程式如下

(3)

式中，Si是i方向(x，y，z)的动量源项，在多孔介质单元中，动量损失对于压力梯度有贡献，压降和流体速度成比例。

湍流模型选择用于描述湍流充分发展的可实现型κ-ε两方程模型，利用基于有限体积法的数值计算方法求解控制方程。研究单相低速流下(即修正的雷诺数Re<10)的流动状态时，黏性阻力的影响较大[9]，在基于多孔介质模型的动量方程源项中设置黏性阻力系数1.7×10-8，内部阻力系数1 300，并按各向异性设置，数量级上相差不超过1 000倍，多孔区域设置为层流模型，其孔隙率保持为1.0。设置入口为速度边界条件，采用均匀流，并通过设置不同速度以求得各速度下的压降值，出口为压力边界条件，采用湍流强度和水力半径方式描述湍流过程，对于湍动能、湍流耗散率等均采用二阶迎风格式的计算方法，多孔介质区域设置为对称(symmetry)边界条件，因为本次物理模型及湍流形式具有镜像对称的特征，在黏性流动计算中也可以使用滑移壁面条件。在进行数值模拟计算过程中，为保证模拟计算良好的收敛性，将流动各个方向速度、湍动能κ和耗散率ε的残差收敛标准均设为0.001[10]。

2.2.5 结果处理分析

结合风洞试验，现从数值模拟软件中选取部分表征压力变化的后处理图(图4、图5)并进行相关分析。

图4 5.5 m/s风速下压力云图

图5 11.2 m/s风速下压力云图

由数值模拟结果和风洞试验数据整理得到入口速度和压降的关系，如图6所示。

图6 压降与速度曲线

2.3 阻沙栅栏数值模拟分析

栅栏作用于地表很大程度上影响了其周围气流的流动特性，是原来流经地表的气流成为一种特殊形式的次生流，进而导致栅栏附近的气流无论是在流动强度还是流动方向都发生了很大的变化，尤其是气流的动量和能量的传递方式变得更加复杂[11]。Plate等最早根据对栅栏防护的空气动力学分析，将绕过二维直立栅栏的流场划分为外层、中间层、内层、中间和外层之间的混合区、内层和中间层的混合区、直立涡旋区和潜在外流区等7个区域[12]。故本文利用基于多孔介质条件的数值模拟软件模拟阻沙栅栏前后流畅的变化，并体现数值模拟软件的适应性和便捷性。

2.3.1 阻沙栅栏几何建模

采用计算机辅助软件建立模型如图7所示，实体栅栏模型如图8所示。

图7 阻沙栅栏几何模型

图8 实体栅栏模型

为了充分描述自然条件下栅栏对来流风速的影响作用，故将计算域取相较于阻沙障碍物尺寸数几十倍大小。此模型沿风速方向取60 m，沿障碍物纵向取0.16 m，计算域高度取10 m，坐标原点位于阻沙栅栏底部中心处。

2.3.2 网格划分

利用CFD数值模拟软件中前处理软件对已有几何模型进行严格的网格划分，最终网格相关数据：558 558个节点、1 623 850个网格面和532 800个网格单元，并在多孔介质区域进行局部加密，网格质量良好，结果如图9所示。

图9 网格划分结果

2.3.3 边界条件和计算参数控制

(1)边界条件的设定

根据空气动力学原理，当马赫数小于0.3时为不可压缩流，故模型入口条件为速度入口，使用UDF编辑入口风速廓线，地表粗糙度设为0.5；对于湍流边界条件，采用湍动能κ和比耗散率ω两方程模型的组合，其中

(4)

(5)

出口和上表面条件因均处在湍流充分发展的情况，故设置为压力出口，湍流计算方法选择为湍流强度和水力半径模式，并设置湍流强度为5，水力半径为0.375；多孔介质区域设置成对称边界条件，并根据以往文献关于压降和风速研究设置相应的分布阻力黏性阻力系数为0.38，惯性阻力系数为300，保持孔隙率为1.0不变。

(2)计算方法及参数控制

湍流计算方法选择定常状态下收敛性较好基于不可压缩流体的SIMPLEC方式，对于动量、湍动能等离散化问题均采用可以减少扩散误差的二阶迎风格式；选择标准初始化从所有区域开始的方法，其他初始值保持缺省状态，为确保模拟结果的准确性，对诸如三个方向的速度、湍动能、比耗散率等的监视，应保证其残差收敛标准控制在0.001，并设空气密度ρ=1.225 kg·m-3，黏度μ=1.789 4×10-5。

2.3.4 计算结果分析

气流在越过栅栏后形成明显减速恢复区，从而降低来流强度，阻止风携沙继续向前运动，障后风速变化最为突出区段主要分布在栅栏后(1～6)H内；现分别取障后z/h=4.2H和y/h=1.38H处风速变化，并与1983年Bradley和Mulhearn等[14]对防风栅栏的风洞试验结果进行对比。见图10、图11。

注：Uo为上风向入口处y=4.0 m高度处的风速 图10 z/h=4.2H风速对比

图11 y/h=1.38H风速对比

由图12分析得知，在栅栏后z/h=4.2H(H为直立栅栏高度)距离处，来流廓线风速度纵向剖面数值与B&M等风洞实验数据整体趋势一致；在y/h=1.38H处，基于多孔介质的数值模拟结果与B&M等在相同条件下相应位置处的风速大小基本一致，相对误差为2.17%。由此可见，利用数值模拟软件中自带的多孔介质模型，通过合理的控制相关参数能有效地模拟一些类似网孔状的实际工程构筑物，减少人力、物力的消耗；再者，从速度云图可以看出，对于非定常来流，阻沙栅栏可以明显减低来流风速，降低来流湍动能，而挤压上扬区并不明显，且在障后形成较不明显的涡流区，减少对障后地表的扰动，此结果与张克存等[15]关于直立栅栏风洞模拟试验或实地观测结果相符。

图12 阻沙栅栏速度云图

直立栅栏作用于地表明显地影响其周围气流的流动特性，气流绕过栅栏经历了抬升、加速、沉降和恢复等过程；栅栏对地表的保护作用表现为风速的削减，这是科研工作者用来评价栅栏防护效益的一个常用指标，当气流越过栅栏后风速急剧降低，最终将导致沙粒减速沉降。本文对直立栅栏基于多孔介质的模拟，从理论上证实了直立栅栏对来流风速的削弱作用和对湍动能的耗散作用，可有效模拟其实际流场变化。

3 结论

采用数值模拟和试验数据结合的方式寻找采用CFD数值模拟软件基于多孔介质模型的最佳边界条件和计算方法，模拟计算中均采用三维建模方法，分别对不同多孔介质模型进行数值模拟分析，结果表明：全断面PE网基于压降的数值模拟，在选择多孔介质模型中，设置动量源项中的黏性阻力系数为1.7×10-8，惯性阻力系数1 300时，可较好地模拟PE网的实际压降变化，且随着入口风速的增加，压降呈明显的指数方式增加。

对阻沙栅栏的模拟表明，多孔介质条件可准确模拟其实际流场变化，通过对计算过程的严格控制，阻沙栅栏可有效降低来流风速，具有显著的防护效果；对于多孔介质模拟的实用性，不同的模型也具有不同的特点，需要根据具体试验或数据分析，确定在模拟过程相关参数的设定，以最合理的方式反映实际情况。

目前多孔介质的应用还没有形成统一理论规范，还望研究者们在模拟的同时多结合风洞试验的结果。

[1] 屈建军，喻文波，秦晓波.HDPE功能性固沙障防风效应试验[J].中国沙漠，2014，34(5):1185-1193.

[2] 康向光，李生宇，王海峰，等.不同组合间距的尼龙阻沙网积沙形态特征对比[J].干旱区地理，2015，38(2):283-291.

[3] 罗万银，董治宝，钱广强，等.栅栏绕流减速效应风洞实验模拟[J]. 中国沙漠，2010，30(1):1-7.

[4] 刘学强，闫晓，肖泽军.多孔介质内单相流阻力特性[J].核动力工程，2009，30(5):40-43.

[5] Albusairi B，Hsu J T. Application of Shape Factor to Determine the Permeability of Perfusive Particles[J]. Chemical Engineering Journal， 2002，89:173-183.

[6] Ergun S. Fluid Flow through Packed Columns[J]. Chemical Engineering Progress， 1952，48:89-94.

[7] Ergun S.Orning A A.Fluid Flow through Randomly Packed Columns and Fluidized Beds[J]. Industrial and Engineering Chemistry， 1949，41:1179-1184.

[8] Teng H， Zhao T S. An Extension of Darcy’s Law to Non-stokes Flow in Porous Media[J]. The Chemical Engineering Science， 2000，55:2727-2735.

[9] 邓彩华，童亮，陈壁峰，等.多孔介质流动的直接数值模拟[J].武汉理工大学学报(交通科学与工程版)，2011，35(6):1257-1260.

[10]赵云鹏，刘兴，董国海.单片网衣周围流场特性的数值模拟[J].中国科技论文在线，2010，5(7):557-562.

[11]罗万银，董治宝，钱广强，等.直立阻沙栅栏流场特征的风洞模拟实验[J].中国沙漠，2009，29(2):200-205.

[12]Plate E J. The aerodynamics of shelterbelts[J]. Agricultural Meteorology， 1971，8:203-222.

[13]陈廷国，马思明.不同截面形式防风网流场的3D数值仿真[J].计算机仿真，2014，31(1):258-263.

[14]E.F.BRADLEY， P.J.MULHEARN. Development of Velocity and Shear Stress Distributions in the Wake of a Porous Shelter Fence. Journal of Wind Engineering and Industrial Aerodynamics， 1983(15):145-156.

[15]张克存，屈建军，牛清河，等.青藏铁路沿线阻沙栅栏防护机理及其效应分析[J].中国沙漠，2011，31(1):16-20.

[16]程建军，庞巧东.戈壁强风区挡风构筑物限制下列车气动力学特性分析[J].铁道标准设计，2013(1):1-5.

Simulation Applicability Analysis of Sand-protecting Barrier along the Railway Based on Porous Media Method

WANG Lian， CHENG Jian-jun， ZHI Ling-yan， XIN Guo-wei

(College of Water Resources and Architectural Engineering， Shihezi University， Shihezi 832003， China)

To simulate the flow field characteristics of multiaperture sand-protecting barrier along the railway with a view to the difficulties in such simulation modeling， the essay uses CFD numerical simulation program to analyze the flow field characteristics and pressure changes of the polythene net and the upright sand-preventing fence based on the 3D simulation of the porous medium condition. With reasonable control of the boundary conditions and parameters， the results show that (1) the numerical simulation of the upright sand-preventing fence is similar to the rate of flow field in wind tunnel experiment， slow recovery area appears obviously after the barrier in porous media simulation because of the effect of rectification of the fence and returns to the same wind speed profile flow at the entrance; (2) combined with the whole section polythene net wind tunnel test， the same simulation results as the test data are obtained by porous media model and the voltage drop increases obviously with the increase of the inlet velocity. The simulation results show that porous media method can better reflect the changes of the flow field of such barrier with proper control of parameters

Porous medium; Sand-protecting barrier; Voltage drop; Flow field

2016-05-10；

2016-05-19

国家自然科学基金项目(51568057；51268050；50908152)

王 连(1992—)，男，硕士研究生，研究方向沙害防治工程，

E-mail:2798895114@qq.com。

1004-2954(2016)11-0018-05

U213.1+54

A

10.13238/j.issn.1004-2954.2016.11.005