被自由面穿过排水孔的数值模拟方法 朱 军,刘光廷,陆述远 摘要:在排水子结构的基础上,结合求解无压渗流场的调整传导矩阵法,给出了被自由面穿过排水孔的数值模拟方法,同时,对渗流量也提出了一种简单的计算方法,最后给出了一个实例分析.排水孔是水利工程中常用的排水设施,它对渗流场水头分布影响很大,如何正确有效地对其渗流行为进行数值模拟就显得十分重要. 但由于排水孔尺寸很小,若将其作为内边界处理,会导致有限元计算网格过于复杂. 目前工程上一般是以给定水头的结点来代替排水孔,这种处理没有反映孔的尺寸效应. 首次提出了排水子结构算法,避免了模拟排水孔时复杂的网格的划分,但无法模拟穿过自由面时排水孔的排水情况,基础上进行了改进,但没有给出排水孔内边界的具体处理方法. 本文在排水子结构的基础上,结合求解无压渗流的调整传导矩阵法,对排水孔与自由面相交的情况给出了具体的算法,有效地解决了排水孔的数值模拟问题,同时对渗流量的计算提出了一种简单的方法. 1 求解无压渗流场的调整传导矩阵法 用调整传导矩阵法求解无压渗流场就是在每次迭代中修正渗透传导矩阵,使其逐步逼近真实的渗流区域. i 次迭代后,对于自由面以上的单元,在下次计算中遗弃;自由面以下的单元,其传导矩阵保持不变;被自由面穿过的单元,其传导矩阵须修正,修正的思路为:单元中水面以上的高斯点渗透系数取原值的千分之一,水面以下的高斯点渗透系数取原值,以此原则重新计算单元传导矩阵.得i + 1 次迭代方程为: [ K]i+1 [ H]i+1 = { Q} i+1 (1) 式中: [ K]i+1 、[ H]i+1 、[ Q]i+1 分别为i + 1 次计算的传导矩阵、待求水头列阵、等效结点流量列阵. 2 排水子结构被自由面穿过时的出口传导矩阵及右端项 对于含有一个排水孔的六面体八结点等参元中的子结构,其内部单元划分(内部分层以2~4 层为宜[ 1 ]) ,如图1. 图1 排水子结构内部单元划分 对子结构可建立如下渗流平衡方程:[ K]z[ H]z = [ Q]z (2) 式中: [ K]z 、[ H]z 、[ Q]z 分别为子结构的传导矩阵、水头列阵、等效结点流量列阵(含子结构邻近单元的结点流量贡献). 子结构传导矩阵由其内部各单元传导矩阵迭加而成,即[ K]z = ∑[ K]ze。 根据上次(设为i 次) 迭代计算出的水头值可判断子结构与自由面的相对位置,当子结构与自由面相交时, 其内部单元高斯点压力水头有的大于零,有的小于零,对于压力水头小于零的高斯点,相应渗透系数取原值的千分之一,压力水头大于等于零的高斯点,渗透系数取原值,以此原则计算各单元传导矩阵,迭加后可得i + 1 次计算的子结构传导矩阵,如式(6) . 当子结构与自由面相交时,其部分内边界(正压力水头内边界,如图2) 排水,可由上次迭代计算出的子结构结点水头值求出排水孔内边界中正压力水头部分的边界流量q ,再由式(5) 等效到排水孔内边界结点上,以修正子结构的等效结点流量列 图2 被自由面穿过的排水子结构阵 次子结构凝聚后的出口传导矩阵及右端项. 自由面以上的子结构,在下次计算中遗弃;自由面以下的子结构,内边界结点水头给定[1 ],传导矩阵及右端项不用修正. 将子结构凝聚后的出口传导矩阵迭加到整体传导矩阵[ K]i+1 上,同时右端项也迭加到整体等效结点流量列阵[ Q]i+1 上,利用方程(1) 计算出结点水头值,由子结构的出口水头值,再利用式(8) 、 (9) 可计算出子结构内部结点水头值,如不满足收敛精度,则按以上思路进行下一轮的迭代计算. 3 计算渗透流量的方法 对于渗透流量的计算,一般方法是先求出计算断面(由于流速在单元面上不连续,计算断面一般由若干单元中截面组成) 上的法向流速,再在断面求积分即得. 本文提出的方法避开法向流速的求解,而仅利用调整传导矩阵法计算的最终传导矩阵与结点水头值加上较简单的运算即可得出渗透流量. 本方法是以若干单元面组成的面为计算断面, 算出计算断面同边的单元(且有一个面在计算断面上) 对断面上各结点的流量贡献,之和即为该断面的渗透流量. 设单元j 的一个面在计算断面上, i 结点为单元j 的一结点且在计算断面上,用下式可计算单元j 对i 结点的结点流量贡献: l Qij = ∑ami, k ×hk (10) k=1 mi 为i结点在j单元中的局部编号, am, k 为j单元最终传导矩阵[ K]je mi 行k列的元素, hk 为j单元水头列阵[ H]je k 行的元素, l 为j 单元的结点个数. 计算断面同边单元对i 结点总的结点流量贡献为:Qi =m ∑Qij , m 为在断面同边与i 结点公点且一个单元 j= 1 面在断面上的单元总数,则通过计算断面的渗流量为: n Q =| Qi| (11) ∑ i = 1 n 为计算断面上的结点总数. 用该方法计算过流量时不用求出法向流速,比常规法简单. 4 算 例 某混凝土重力坝,坝高100 m , 底宽70 m , 顶宽10 m , 上游水位为89 m , 下游水位为11 m , 在坝体中距上游面3 m 处设孔径为15 cm 、孔径为5 m 的垂直排水孔幕,坝基也设有深30 m 的垂直排水孔幕,其孔矩、孔径与坝体排水孔相同,主要几何尺寸如图3. 图中A -A、B -B 为计算渗流量的两个断面. 坝体、坝基的渗透系数分别为:1. 0 ×10-8 、1. 0 ×10-7m/s. 有限元计算网格如图4 , 计算区域沿厚度方向取5 m , 坝体、坝基各设一排水孔,排水孔的中心在厚度方向2. 5 m 处. 图3 计算区域主要尺寸 图5、6 为不设排水孔时水头等值线图,图7、8 为设排水孔时沿排水孔中心断面水头等值线图. 由图6、8 对比可知,在坝体,由于排水孔的作用,使自由面位置大幅下降,并使坝体中的渗流水头主要集中消耗在排水孔附近,降压作用十分显著. 由表1 可知排水孔的排水作用也是很明显的,在不设排水孔时,A -A 断面的渗流量为:3. 478 998 ×10-5m3/ s , 而设排水孔后其渗流量减小为2. 028 225 ×10-5 m3/ s , 其余渗透水量通过排水孔排出了坝外. 由图9 可以看出坝基排水孔的设置有效地降低了坝基面的扬压力,这对维持坝体整体稳定十分有利.图5 无排水孔时坝体坝基等水头线 图4 有限元计算网格图 6 无排水孔时坝体等水头线图7 有排水孔时排水孔中心断面坝体坝基等水头线 图9 坝基有无排水孔时坝基面水头分布对比(沿排水孔中心断面) 5 结 语本文对排水孔与自由面相交的情况给出了具体的算法,有效地解决了排水孔的数值模拟问题, 实例的计算结果符合工程实际规律. 文中同时对渗流量的计算也提出了一种简单的方法,该方法的计算结果满足渗流连续性要求. 这些研究成果对工程有一定的实用价值. |