发布时间 : 星期二 文章用动网格模拟闸门开启过程非恒定水流特性更新完毕开始阅读
用动网格模拟闸门开启过程非恒定水流特性 沙海飞1，2 周辉1，2 吴时强1 陈惠玲1
（1.南京水利科学研究院 南京 210029；
2.水利部水科学与水工程重点实验室 南京 210029）
【摘 要】 本文采用非结构动网格技术，采用网格变形与局部重构相结合的方法实现动网格，网格的变形采用弹簧模型，对闸门开启过程非恒定水流特性进行了数值模拟，得到了不同开启速度时的流量、压力分布以及不同流态的特征与发生、发展过程，结果令人满意。 【关键字】 非结构动网格 闸门 非恒定水流 数值模拟
0 前言
在水工建筑物中，闸门是用于控制水流最为常用的设施。在水坝、船闸等设计过程中，工程技术人员必须了解上下游水流的运动规律。由于闸门前后存在水位差，当水坝、船闸
[1]
的闸门开启时，门后容易发生水流空化，时常造成空蚀破坏与闸门振动等危害。文献[2]发现灌区水闸在闸门开启过程中，过闸水流为非恒定流，闸门开启的时段内闸下产生的冲坑深度，占一次放水过程产生冲坑深度的一半。可见在闸门开启过程非恒定水流是要引起足够重视的。而目前这方面的研究主要结合具体工程进行缩尺模型试验，闸门开启时的非恒定水流与复杂的漩涡运动增加了试验的难度，一般通过测量壁面压力来推断体型是否合理，对闸门后的内部流态了解很少。
随着计算机技术的飞速提高，计算流体力学得到了迅速发展，采用数值模拟闸门开启过程非恒定水流特性成为可能。文献[3]对假设闸门突然开启的情况进行了数值模拟，这和实际偏差较大。文献[1]对船闸弧门非恒定流态进行了二维计算，其计算网格是固定的，对弧门处网格进行了特殊处理，和实际也有一定的偏差。文献[4]用动网格技术对挑流冲刷过程进行了数值模拟，采用的是结构网格，用网格变形来实现动网格，这种方法很难用于闸门开启这种复杂的运动过程。随着非结构动网格技术的发展，以及在军事、仿生等领域的成功应用，使精确模拟闸门开启过程非恒定水流的特性成为可能。
本文采用非结构动网格技术，对平板闸门不同开启速度时非恒定水流进行了垂向二维数值模拟，探讨了闸门后流态发生、发展机理。
1 动网格技术
本文所处理的问题位移较大，所以采用网格变形与局部重构相结合的方法实现动网格。基本做法是，开设一个包围运动边界的窗口，窗口内的网格可以变形以适应边界的运动。以单元面积为判据，当窗口内即变形区出现严重扭曲的网格单元时，删除变形区内的网格单元，然后用阵面推进法重新生成该区网格。流动参数通过线性插值由旧网格映射到新网
中国水利学会第二首届青年科技论坛论文集
312
格。
控制网格变形的弹簧近似模型是将变形区网格看作一个弹簧网络系统，每一条边都认为是一根具有一定倔强系数的弹簧[5]，网格点所受的合力为 fi??K?xijj?1Nij?xi? (1)
其中Kij为连结节点i、j的弹簧的倔强系数，xi是节点i的位置矢量，Ni是与节点i相连的节点数。网格点的受力始终等于初始状态所受的合力，上式可表示节点i的初始受力，则当运动边界上的点移动后，变形区内节点的新位置坐标可以通过迭代求解如下的线性系统得到。
???K1j?j?? ???????K2jjaiki?k???x??f?11??x??f?aik?i?k?2???2? (2) ??????????xf?sumsum??????Ksumj?j??0，若k点不与i点相连 aik?? (3)
K，若k点与i点相连?ik该系统对角占优，Jacobi迭代格式为： xk?1i?Ni????Kijxk?fji?j?1???Kj?1Niij (4)
当边界运动时，下个时刻的固定边界点和主动边界点的位置解已知，因此弹簧系统的边界
条件为Dirichlet型，上式经过数次迭代就可达到满意的精度。
2 控制方程与数值求解方法
考虑是不可压缩水流流动问题，基本控制方程为Reynolds方程，省略各变量的时间平均项标志“”，紊流模型采用k??双方程模型，方程如下：
连续方程： ?ui?0 (5)
?xi? (6) 动量方程： Dui?p???uiDt?f????uiuj???x?i?x?x?j?ij?k方程: Dk????????k?????t???Gk?? (7)
Dt?xi???k???xi??????t??????2 (8) 方程： D????????C1?pk?C2??????x?Dt?xi??kk??i????式中，t为时间；ui和xi分别为速度分量和坐标分量；?为运动粘性系数；p为修正压力；
fi为质量力；?t＝cuk2?为紊动粘性系数；方程中的经验常数cu?0.09，?k?1.0，
???1.33，C1??1.44，C2??1.42 Gk为平均速度梯度引起的紊动能产生项，
313 四、水力学 河流动力学
??uGk＝?t?i?x??j??ui?uj??；对于雷诺应力项?uu的封闭问题，这里用k??双方程模型来处
ij??????xj?xi?????理。
采用有限体积法对方程进行离散，对时间和空间均采用二阶精度格式，空间离散采用中心差分格式，对时间采用二阶精度的离散格式是很又必要的，因为研究的闸门开启过程水流随时间的演变是很重要的。压力校正采用SIMPLEC算法，并采用显式校正步法对SIMPLEC算法进行了显式校正。
上下游均是压力边界条件，给定上下游的总压，对于很长管段中闸门区域的局部模拟，可以采用一维计算为二维计算提供边界条件；壁面采用无滑移边界条件，采用壁面函数来处理壁面边界。
3 实例计算及分析
本文计算区域如图1所示，模拟闸门开启过程非恒定水流特性，上下游水深分别为30 m和10 m，闸门长6 m，开启速度为v(假设闸门是匀速开启的)，点A和B分别是闸后1.5 m处上下壁面处的压力测点，计算中不考虑上下游水位的变化。计算时间步长为0.01s，计算区域的网格每个时间步长调整一次，经过单个时间步长的网格位移比较小，采用网格变形的方法处理，经过数步之后，网格扭曲就比较厉害，进行局部重构，图2是开启时间为10s的不同时刻的网格情况，网格质量基本较好。共模拟了开启时间为10s、20s、30s和60s四种情况，也就是开启速度分别为0.6m/s、0.3m/s、0.2m/s和0.1m/s。 30640B0.510010Av t=6.67s t=8.33s t=10.00s 图1 计算区域示意图 t=0.00s t=1.67s t=3.33s t=5.00s 图2 不同时刻网格（开启总时间为10s） 3.1 流量变化规律 图3是不同开启时间的开启过程流量曲线图，在闸门开启前，其流量都为零，随着闸门的逐渐开启，流量都随之增大，在达到全开之后，流量都有一个平稳的增长过程。不同开启时间在开启过程中流量的变化有一定的差异，开启时间短的过程基本曾线性增加，开 中国水利学会第二首届青年科技论坛论文集
314
启时间变长后流量的波动明显加剧，这是由于开启速度的减缓，形成了较大尺度的漩涡，对流量形成了明显的影响。不同开启时间在达到全开时的流量也有明显的差异，开启时间短的过程在达到全开时的流量较小，之后流量的增加也较快，而开启时间长的过程在达到全开时的流量较大，之后的增速也较缓，慢慢接近稳态时的值。 120100Q(m3/s·m)10s 20s 30s 60s 80604020001020304050607080t(s) 图3 开启过程流量曲线图 3.2 流态特征分析
图4是开启时间为10s、20s和30s的不同时刻流态图，为了便于比较，3种开启时间列出的相对开度是一致的。不同开启时间在较小开度时，闸门后回流区随开度逐渐增大而逐渐形成、发展，在开度增加的过程中，主流厚度在增加，并不断的将原先处于主流上侧的低速水体加速汇入主流，从而能在启门后一段时间内维持门后漩涡的稳定，开启速度较慢时形成的漩涡较大较稳定。随着闸门的进一步开启，闸门后漩涡出现分裂、脱落，接着
Frame001?08Aug2005?shahaifei_10s新的漩涡生成、运动、合并、分裂、脱落等引起水流较剧烈变动，在门后不断有尺寸较小Frame001?08Aug2005?shahaifei_10sFrame002?08Aug2005?shahaifei_20sFrame003?08Aug2005?shahaifei_30sFrame002?08Aug2005?shahaifei_20sme001?08Aug2005?shahaifei_10sFrame003?08Aug2005?shahaifei_30s的次生涡产生，随着开度增大，门后回流区尺寸随开度增加而逐渐减小，漩涡产生分裂的Frame002?08Aug2005?shahaifei_20s10Frame003?08Aug2005?shahaifei_30st=1.0s,e=0.101010周期也不断缩短。 t=1.0s,e=0.10t=2.0s,e=0.1010101010100t=1.0s,e=0.1010t=2.0s,e=0.100t=3.0s,e=0.105t=1.5s,e=0.15t=3.0s,e=0.15t=2.0s,e=0.100t=3.0s,e=0.105t=3.0s,e=0.100500t=1.5s,e=0.15t=3.0s,e=0.150t=1.5s,e=0.15t=3.0s,e=0.15t=4.5s,e=0.15t=2.2s,e=0.22t=4.4s,e=0.22t=6.6s,e=0.22t=4.5s,e=0.45t=9.0s,e=0.45t=13.5s,e=0.45t=6.5s,e=0.65t=13.0s,e=0.65t=19.5s,e=0.65t=7.0s,e=0.70t=14.0s,e=0.70t=21.0s,e=0.7000-10-5-10-10-5t=4.5s,e=0.15-10t=2.2s,e=0.22t=4.4s,e=0.22t=4.5s,-10e=0.15-5t=2.2s,e=0.22-10t=4.4s,e=0.22-10-10-10-20-15-20-20-25-30-30-20t=6.6s,e=0.22-15t=4.5s,e=0.45-20-20t=6.6s,e=0.22-15t=4.5s,e=0.45-20t=9.0s,e=0.45-20t=13.5s,e=0.45-30t=6.5s,e=0.65-25t=13.0s,e=0.65YYYt=9.0s,e=0.45YY-20Y-30t=13.5s,-30e=0.45t=6.5s,e=0.65-25t=13.0s,e=0.65YYY-30-30t=19.5s,-40e=0.65t=7.0s,e=0.70t=14.0s,e=0.70-30t=19.5s-40,e=0.65t=7.0s,e=0.70-30-40-35-40-40 -35-40t=14.0s,e=0.70t=21.0s,e=0.70-35-4020-4020-456020304050-40-45202030XX(a) 开启时间 70 (b) 20 30 80(c) 30s -45405010s 60802040开启时间20s 6040开启时间506040602020t=21.0s,e=0.70 30 403050 70 80 4070X60708060 70808040506080图4 闸后流态分布图 803040506070803.3 压力分布特征 在闸门开启过程中闸门前后的压力变化较大，图5是不同开启时间情况下，闸后A点和B点的压力随相对开度变化的情况。在e?0时A点和B点压力为静水压，压力分别为10×9.8kPa和4×9.8kPa。在闸门开启的瞬间，压力陡然增大，开启速度越快压力增大值越