作者：郑 俊 杨晓娟 王昱倩 金鑫鑫 费文平
论文关键词：无压；渗流；自由面；数值计算
论文摘要：在水利水电工程中，存在许多有自由面的无压渗流问题，自由面是渗流场特有的一个待定边界，这使得应用有限元法求解渗流场问题时，较之求解温度场和结构应力等问题更为复杂。归纳总结了无压渗流分析的各种数值计算方法，分析比较了其优缺点和适用条件，提出了无压渗流数值分析方法的发展趋势。
1 引言
在许多水利工程中(如土石坝渗流、混凝土坝渗流、拱坝绕流、地下结构渗流等等)，都存在着无压渗流问题，这类问题的关键在于求解渗流场的边界，即确定事先不知道其位置的自由面和溢出面，属于非线性边界问题。求解该问题的有限元法以往采用移动网格法。虽然取得了许多成功的经验，但也表现出方法本身的缺陷。为解决上述问题，国内外学者致力于寻找有自由面渗流分析的新方法。其研究核心就是计算中不变网格，自neumann于1973年提出用不变网格分析有自由面渗流的galerkin法以来，出现了多种固定网格法，如剩余流量法、单元渗透矩阵调整法、初流量法、虚单元法和虚节点法等。WWw.11665.Com
2 无压渗流的数值分析方法
2.1 调整网格法
调整网格法先根据经验假定渗流自由面的位置，然后把它作为一个计算边界，按照vn=0的边界条件进行分析，得出各结点水头h值后，再校核h=z是否已满足。如不满足，调整自由面和渗出点的位置，一般可令自由面的新坐标z等于刚才求出的h，然后再求解。
该方法原理简单，渗流自由面可以随着求解渗流场的迭代过程逐步稳定而自行形成，并且迭代是收敛的。但是，当初始自由面与最终自由面相差较大时，容易造成迭代中的网格畸形，甚至交错重叠；当渗流区内介质的渗流系数不均匀时，特别是有水平分层介质时，程序处理困难；对复杂结构问题，由计算机自动识别和执行网格移动几乎是不现实的。
2.2 剩余流量[1]
剩余流量法通过不断求解流过自由面的法向流量（称为剩余流量）建立求解水头增量的线性代数方程组，达到修正全场水头和调整新的自由面位置的目的。迭代过程中只需一次形成总体渗透矩阵，但需要判断自由面被单元分割的各种情形，要求算出穿过单元的自由面被单元切割的面积及流过自由面的法向流速，计算工作量很大，难以推广到三维问题中。剩余流量法的全部调整均基于第一次有限元计算的结果，因而计算精度较差。
2.3 单元渗透矩阵调整法[2]
单元渗透矩阵调整法利用对渗流场有限元计算的结果，根据单元结点水头与结点位置势的比较，把渗流场进行分区，各区的渗透系数给不同的值，通过不断调整单元渗透矩阵，模拟渗流不饱和区的作用，来确定出真实的渗流饱和区及渗流场。
该算法实际上是把边界不确定的非线性问题转化成了材料非线性问题来考虑。但是，单元渗透矩阵调整法对三维而言其计算效率是很低的，不能真实反映渗透区域的透水特性，计算精度和收敛稳定性都受到影响。
2.4 初流量法[3]
初流量法利用高斯点的水头求出结点的初流量作为求解水头增量的右端项，避免了求自由面被切割的面积，同时避免了每次迭代中确定自由面的位置的做法，大大简化了剩余流量法的计算工作量。由于初流量法在计算跨自由面单元的结点初流量时，自由面以下的高斯点未予计算，计算精度受到影响。初流量法其收敛性不尽人意，解的稳定性不好。
2.5 虚单元法[4]
虚单元法以上一次有限元计算的结点水头值为基础，求出自由面与单元边线的交点，移动跨自由面单元的某些结点，使之落于交点处，自由面将单元分成渗流实区和虚区。渗流虚区在下一次计算中退出计算区域，随着渗流计算区域向渗流实区逼近，结果也逼近问题的真解。该方法对三维复杂问题不适用，易产生结果收敛不稳定的现象。同时，虚单元法在处理有自由面穿越的单元时，结点移动路径的确定是比较困难的。
2.6 虚节点法[5]
虚节点法以上一次有限元分析求得的节点势为基础，求出自由面和单元节线的交点，根据交点确定单元的积分区域，形成下一次分析的渗透矩阵。不同于虚单元法，虚节点法无需移动任何节点，因此不会出现网格畸形；虚节点法对网格不作改动，并能精确地描述跨越自由面单元的渗透矩阵，具有很好的精度和数值稳定性。
此外，无压渗流的数值分析方法还有边界单元法、流形单元法、无单元法等。
3 无压渗流数值分析方法的比较
调整网格法计算原理简单，迭代过程稳定而自行形成，迭代过程收敛，但该算法对有复杂夹层和复杂排水系统的水工结构处理起来太困难，几乎不可能实现；另外对初始渗流自由面位置的假定要求也较高，如果初始位置与最终自由面位置相距甚远，则极易造成单元严重畸变，影响计算的精度；剩余流量法计算工作量很大，难以推广到三维问题中。初流量法在剩余流量法的基础上作了重大改进，大大简化了剩余流量法的计算工作量，但是收敛稳定性较差，而且由于两种算法的整个迭代过程依赖于第一次有限元计算的结果，精度受到一定的影响。单元渗透矩阵调整法对跨自由面单元按复合材料单元处理，复合材料单元渗透系数在复合面突变，其单元渗透矩阵不能代表这一特性，且矩阵主系数常不占优，因而计算精度和计算稳定性均受到影响。虚单元法对三维复杂问题不适用，易产生结果收敛不稳定的现象。虚节点法具有很好的精度和数值稳定性。
结论
本文归纳总结了各种无压渗流数值计算方法的原理及其优缺点，得到如下结论：
传统的调整网格法虽仍被使用，但由于自身的缺陷给应用带来诸多不便，因而正在逐渐被固定网格法所取代。具体选择计算方法时，应从问题的复杂度、收敛性及精度要求等方面加以考虑。现有的大型商用软件如ansys提供了良好的二次开发环境，用户可以通过二次开发，来实现无压渗流的数值分析。
参考文献
[1] desai c s. finite element residual schemes for unconfirmed flow [j]. int num method eng. 1976, 10(6):1415~1418.
[2] bathe j n. transmit matrix method for seepage with free surface problem [j]. int j num meth engng, 1983, (7):41~53.
[3] 张有天, 陈平, 王镭. 有自由面渗流分析的初流量法[j]. 水利学报, 1988, (8):18~26.
[4] 吴梦喜, 张学勤. 有自由面渗流分析的虚单元法[j]. 水利学报，1994, (8):67~71.
[5] 凌道盛. 有自由面渗流分析的虚节点法[j]. 浙江大学学报(工学版), 2002,36(3):244~246.
