利用Fluent软件对植物坝覆盖区域流场的数值模拟

刘 锋,邱秀云,周 著

(1.新疆农业大学发展规划处,新疆乌鲁木齐 830052;2.新疆农业大学水利与土木工程学院,新疆乌鲁木齐 830052)

1995年,毕慈芬、李桂芬和于倬德针对我国砒砂岩地区水土流失对沟道径流产粗沙以及暴雨洪水侵蚀单一的特点,提出采用沙棘植物“柔性坝”来拦截沟壑粗沙[1]。历经数年,已将其部分粗沙固于沟壑,减少了进入黄河的粗沙,同时也使这些地区的生态环境得以改善。所谓沙棘植物“柔性坝”(以下简称植物坝)就是在沟道或河渠主槽内,按一定的株距和行距种植一定数量的沙棘植物。经调查,砒砂岩在新疆普遍存在,而且因其生存环境恶劣,产生的粗沙相对较多,因此更迫切需要解决。

植物坝的研究为西部地区,尤其为新疆地区的防沙、治沙和水土流失治理工作提供了新思路。新疆农业大学水利与土木工程学院水工实验室于2002年和2003年进行了植物坝研究的室内水槽清水模型试验,分别探讨了植物坝在行列对齐及交错布置方式下,不同种植密度对河渠水流的影响[1-2]。为了进一步研究植物坝的固沙机理,笔者所在的课题组进行了室内水槽浑水模型试验,并利用Fluent软件对植物坝覆盖区域的流场进行了数值模拟。

1 物理试验简介

1.1 试验设备

试验在新疆农业大学水利与土木工程学院水工实验室的活动玻璃水槽内进行(图1)。水槽长22m,宽0.4m,高0.5m,槽首设有进水阀门和静水池,静水池中设消能栅板,用于稳定水流。水槽末端设有可调节尾门,水流经尾门后跌入水箱。水箱出口为一矩形薄壁堰,用于量测流量。水槽上设有活动测针架,可测读水位、床面高程等。水槽顶部有带刻度的轨道,可供测针架移动和读取所测断面位置。将水槽轨道零点刻度所在的断面记为SC0+00断面,1m刻度所在的断面记为SC1+00,其余断面根据所处位置的刻度依此类推。水槽中端底部设一固定铰,首、尾端底部下分别设活动铰各1座,用于调节水槽底坡。

试验中,水位测量采用精度为0.1mm的水位测针量测,流速测定选用南京水利科学研究院生产的LGY-Ⅲ型多功能智能流速仪。

图1 水槽布置示意图

1.2 模型比尺

该研究的主要目的是为了解决新疆河流的多粗沙问题。通过野外调查,并且考虑便于研究,主要依据模型树的高度来确定长度比尺λL=10,选用1∶10的正态模型,采用中值粒径d50=0.34mm的非均匀天然沙为推移质模型沙。

1.3 模型树的选择

为了较真实地模拟天然沙棘对水沙的影响,同时考虑树木的阻力和变形相似,该试验采用装饰用的绿色塑料树作为模型树。塑料树全高约10cm,主干直径约2mm,上部有树枝和树叶,树顶最大横向宽度为7cm,并具有一定的柔韧性。模型树的根部有可用于固定的圆形孔,试验时将模型树固定在有插座的网格上。模型树排成清水试验中优选出的阻水效果最好的6cm×9cm×9cm行列交错布置方式,如图2所示,种植棵数为63棵,种植密度为248棵/m2。

图2 模型树布置方式示意图(单位:cm)

2 流场的数值模拟

水流在植物坝内部的流动可以看作是一种多孔介质内部的流动[3]。无论是从物理的角度还是从数学的角度看,多孔介质内部的孔隙流都是一种十分复杂的流动现象,本文做出如下假定:①多孔介质为各向同性、均匀的刚性介质;②多孔介质上部边界为刚性可渗边界。

很多问题中包含多孔介质的计算,比如流场中包括过滤纸、分流器、多孔板等边界时就需要使用多孔介质条件。在计算中可以定义某个区域或边界为多孔介质,并通过参数输入,定义通过多孔介质后流体的压力降。

2.1 多孔介质模型的假设和限制条件

多孔介质模型采用经验公式定义多孔介质上的流动阻力。从本质上说,多孔介质模型就是在动量方程中增加了一个代表动量消耗的源项。因此,多孔介质模型需要满足下面的限制条件:①因为多孔介质的体积在模型中没有体现,在缺省情况下,Fluent软件在多孔介质内部使用基于体积流量的名义速度来保证速度矢量在通过多孔介质时的连续性。如果希望更精确地进行计算,也可以让Fluent软件在多孔介质内部使用真实速度。②多孔介质对湍流的影响仅仅是近似。③在移动坐标系中使用多孔介质模型时,应该使用相对坐标系,而不是绝对坐标系,以保证获得正确的源项解。

2.2 多孔介质模型对湍流的处理

Fluent软件在多孔介质计算中通过求解标准守恒型方程计算湍流变量。在计算过程中,通常假设固体介质对湍流的生成和耗散没有影响。在多孔介质渗透率很大,而介质的几何尺度对湍流涡结构没有影响时,这个假设是合理的。在另外一些情况中,可能不需要考虑流动中的湍流,即假定流动为层流。如果计算中使用的湍流模型是或Spalart-Allmaras模型中的一种,则可以通过将湍流黏度 μ设为零的方式忽略湍流的影响。如果 μ被设置为零,则在计算过程中仍然会将湍流变量输运到介质的另一面,但是湍流对动量输运过程的影响则完全被消除。在流体面板中将多孔介质区设为层流区选项就可以将 μ设为零。

有限厚度的多孔介质压强变化是用Darcy定律和一个附加的惯性损失项来定义的:式中:μ为湍流黏度;a为渗透率;ρ为水的密度;ν为法向速度;C2为压强跃升系数;Δm为介质厚度。

多孔跃升模型的相关参数在多孔跃升面板中输入。需要输入的参数如下:多孔跃升区域,a,Δm,C2。如果需要考虑弥散相,在多孔跃升区域定义弥散相边界条件。

本文在Gambit里将多孔介质设置为fluid,在fluent里设置多孔介质参数,压强跃升系数为零,这样可以保证流过多孔介质的层流中,压力降正比于速度,同时可以通过试算来确定a值。

2.3 网格剖分

网格生成是计算流体力学的一个重要分支,数值计算中用离散的网格来代替原物理问题中的连续空间,网格中的节点则是所求解物理量的几何位置[4]。需要先对计算区域进行网格划分,而后在划分的网格上对控制方程进行离散、求解。

网格质量对计算精度和稳定性有很大的影响。网格质量包括节点分布、光滑性以及网格变形程度。

本文采用较常用的Gambit软件来生成网格,选用结构网格对模型进行分区网格划分,这样保证了网格质量,减少了整体网格单元数目,同时可以较好地验证物理试验中以某个断面为对象进行测量的结果。植物坝覆盖区段共生成网格单元约2800个,网格节点约2911个。其他区域网格单元约4 641个,网格节点约7200个。

2.4 边界条件的设置

在数值模拟计算域内的流动是由边界条件驱动的,从某种意义上说,求解实际问题的过程,就是将边界线或边界面上的数据根据控制方程外推扩展到计算域内部的过程。因此,提供符合物理实际且适当的边界条件极其重要。边界条件设置时必须保证在合适的位置选择合适的边界条件,同时让边界条件不要过约束,也不要欠约束。本文水气两相流场计算所涉及的边界条件有:进口边界条件、出口边界条件、固壁边界条件和自由水面。

2.4.1进口边界条件

进口边界包括:进水口边界(图3);模型顶部为空气进口边界。

图3 边界条件示意图(纵向剖面)

进水口全部被水充满,可以看做是均匀来流,采用速度进口边界条件,在2个物理模型中均以进口水流的方向为 x轴正向,以竖直向上为y轴正向。边界条件表达式为式中:ux,uy,uz分别为x,y,z方向的分速度;U0为来流速度,根据实测入流量而定。

空气进口边界采用压力进口边界条件,压强设定为大气压,即p=0。全部为气体。

2.4.2出口边界条件

出口边界为压力出口和自由出流。出口由水和空气两部分组成,压力出口压强设定为大气压,即p=0;自由出流完全发展的区域:

式中:n为出口截面外法线方向;k为紊动能;ε为紊动耗散率。

2.4.3固壁边界条件

在固壁边界上采用无滑移边界条件,即并利用实测的数据确定固壁糙率。

2.4.4自由水面

采用VOF方法捕捉水气交界面,得到自由水面。

对于边界上的紊流参数k和ε,其初始值如果不能从试验中取得,只能采用一定的经验公式来估算。由式(5)(6)来估计进口的k和ε。式中:I为底坡;u0为入口处的速度;Cμ为跃升系数,可取0.09;l为湍流长度;L为关联尺寸,对于充分发展的湍流,可取L等于水力半径。

2.5 数值计算过程

本文模型的数值计算过程分2步,先是模拟水槽内无植物坝覆盖时水流定常流流动,直至收敛,然后在这个收敛流场基础上模拟水槽内有植物坝覆盖时的水气两相非定常流流动。否则如果直接模拟水气两相非定常流,迭代计算容易发散[5]。在此基础上应用VOF方法捕捉自由水面,变成水气两相非定常流模拟。数值模型的边界条件设置:水流入口采用速度入口,水流入口水位根据试验时的实际工况设定;入口水位以上至下游出口断面设置为气体压力进口,静压设为零;下游出口分为空气出口和液态水流出口两部分,空气出口部分采用压力出口,静压为零;各边界条件的紊流参数均按照经验公式确定湍动能和紊动耗散速率,固体边壁当量粗糙度以试验取0.05m。压力插值格式选用body force weighted格式,动量、紊动能、紊动耗散速率的离散均采用一阶迎风格式[6]。时间步长为0.001s,需求时间步设定为100000,每个时间步内的最大迭代计算次数为20次。

图4为模拟有植物坝覆盖水槽的3个典型断面上的流速值与试验值的对比。

图4 模拟流速值与试验值的对比

水流通过植物坝前端时,随着植物坝阻水所产生的壅高水流流速急剧减小;水流流至植物坝中段时,流速趋于稳定,沿水深纵剖面上流速值变化不大;水流通过植物坝末端时,流速受植物枝叶对水流反作用力的影响,随着水深的增加流速值再次减小。

由图4中3组典型断面数值模型计算流速值与物理模型试验值的比较可以看出,3个断面上的模拟流速沿水深分布趋势与试验测得的结果基本一致。由此表明,本文利用Fluent软件计算和VOF方法耦合的数学模型是正确的,用多孔介质域模拟植物坝覆盖区域是可行的。

3 结 语

a.水槽内放置植物坝后,槽内水流的流场与植树前有了明显的不同,其中上游壅水区是植物坝对上游水流阻挡作用的直接结果。

b.将数值计算得到的流速值与模型试验流速值实测结果进行了比对,两者吻合较好。表明本文的数学模型是正确的,用多孔介质域模拟植物坝覆盖区域是可行的。

由数值模型计算流速值与物理模型试验值的比较表明,Fluent软件计算和VOF方法耦合的数学模型是正确的,用多孔介质域模拟植物坝覆盖区域是可行的[7]。通过较系统的有植物坝河道的试验研究分析,得到了一些有价值的数据和结论,但这仅仅是初步的研究工作,还有待于深入的理论分析和采用更多的数学工具进行更多的物理模型试验和数值模拟方面的研究。

[1]程艳,李森,邱秀云,等.河渠种树水流特性试验研究[J].新疆农业大学学报,2003,26(2):59-64.

[2]阎洁,周著,邱秀云,等.植物“柔性坝”阻水试验及固沙机理分析[J].新疆农业大学学报,2004,27(3):40-45.

[3]袁梦,黄本胜,邱秀云,等.水葫芦覆盖区水流阻力效应试验研究[J].广东水利水电.2008(2):29-34.

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

[5]陈界仁,朱晓丹.含柔性沉水植物明渠水流运动特性试验研究[J].水利水电科技进展,2008,28(3):16-19.

[6]吕升齐,唐洪武,闫静.有无植物条件下明渠水流紊动特性对比[J].水利水电科技进展,2007,27(6):57-60.

[7]吴航,拾兵,陈志娟.淹没状态下河道植物偏斜高度与糙率的关系[J].水利水电科技进展,2008,28(1):20-22.