河道平面二维水沙数学模型的有限元方法.docVIP

　　河道平面二维水沙数学模型的有限元方法 摘要：采用有限元 方法 建立起一套河道平面二维水流泥 沙数学模型。在前人 研究 的基础上，采用了质量集中的处理方法，提出了压缩存储的方法， 从而大大减少了 计算 存储量。针对有限元法时间步长需取得较短 问题 ，采用了“预报-校正-迭代”的算法，提出了“非恒定-恒定-非恒定流”的算法，既能解决工程实际问题，又大大减少了计算量。作者以下荆江监利河段为例进行泥沙冲淤计算，计算结果与实测值符合较好，从而证明了模型的可靠性。 关键词：水流泥沙 有限元 模型验证 三峡工程建成后，水库将拦蓄大量泥沙，下泄水流含沙量减小，对三峡工程坝下游河道将产生以冲刷为主的 影响 ，包括对荆江河段的河势及荆江大堤带来影响。为研究坝下游重点河段的河床冲淤分布、河势变化、近岸流速变化等问题，一维模型显得无能为力，但可采用平面二维模型来解决。有限元方法可采用无结构化网格，能很好地模拟不规则的几何形状，因此很适合于对天然河道的模拟。 然而，正如其它方法一样，有限元法也有它的缺点，主要是计算存储量和运算量较大。为扬长避短，使有限元方法能运用到对天然河道的模拟上来，本模型运用质量集中［4］的方法将系数矩阵转化为三对角矩阵，并提出了紧凑的分块压缩存储方法，从而大大减少了计算存储量，使得计算能在一般微机上进行。采用质量集中方法的不足之处是时间步长需取得较短，且在河道模拟中尤为突出(因河道比较窄长，网格需划分很细，而该法的稳定性要 求时间步长与网格尺度成正比)。针对该问题，笔者采用了“预报-校正-迭代［5］”的算法，该法可加大时间步长，同时有效避免了数值震荡。针对长系列水沙条件下计算量较大问题，作者又提出了“非恒定-恒定-非恒定流”的算法，该算法既能解决工程实际问 题，又大大减少了计算量，使有限元方法能够很好地运用于河道水流泥沙问题的实际计算。 1 基本方程 　　平面二维水流方程 　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　(1) 　　　　　　　　　　　　　　　　　　　　(2) 　　　　　　　　　　　　　　　　　　　　(3) 悬移质泥沙扩散方程  　　　　　　　　　　　　　　　　　　　　(4) 推移质不平衡输移方程［6］ 　　　　　　　　　　　　　　　　　　　　　　　　　(5) 河床变形方程 由悬移质引起的河床变形方程为 lt;![endif]gt;　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　(6) 　　由推移质引起的河床变形方程为 　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　(7) 以上各式中U，V分别为垂线平均流速在x,y方向上的分量；Zs、Zb和H分别为水位、河底高程和水深；g为重力加速度；vt为水流紊动粘性系数；ρ为水的密度；τx、τy、舄瓂分别为底部切应力在x和y方(τx、τy)=，向上的分量：C为谢才系数，常用曼宁公式计算：C=H1/6/n；S和S*分别为垂线平均含沙量和挟沙力；N和N*分别为推移质输沙量和推移质输沙能力折算成全水深的泥沙浓度；εs为泥沙紊动扩散系数；ω为泥沙沉速；γ′为床沙干容重；α为悬移质泥沙恢复饱和系数，淤积时取0.25，冲刷时取1.0；β为推移质泥沙恢复饱和系数，取0.25。 悬移质泥沙共分成8组，水流挟沙力和分组挟沙力级配采用李义天方法［7］进行计算。 2　有限元方程 将整个计算域剖分成一个三角形网格系统。每个三角形为一个单元，其编号为e,e=1,2,3， …，NE,NE为单元总数。单元三个顶点为节点，其局部编号为j=1,2,3(以逆时针为序)。节点的整体编号为i,i=1,2,3,…，NP，NP为节点总数。节点的整体编号与局部编号计算前一定要规定好。 引入插值函数：f=fi(t)φi(x,y),φ为形函数。对方程(1)～(5)中的变量用插值函数近似表示，并使用伽辽金法［1］对方程进行整理变形，可得到积分方程 　　　　　　　　　　　　　　　　(8) 　　　　(9) 　　　　　（10） 　　　　　　(11) 　　　　　　　　　　　　　　(12) 由此，经整理得到如下有限元方程 　　　　　　　　　　　　　　　　　　　(13) 　　　　　　　　　　　(14) 　　　　　　　　　　　　(15) 　　　　　　　　　　　(16) 　　　　　　　　　　　　　　　　(17) 　其中 i,j,k=1,2,3,…，NP. 3　数值解法 　　为书写方便，采用通用变量P来代替方程(13)～(17)中的变量Zs、U、V、S和N，用FP表示方程中等号右边项，并用对角矩阵 来代替Aij，则方程可化成统一的形式： dPi/dt=FPi。 对该方程的求解，模型采用了“预报-校正-迭代［5］”的 计算 方法 。其方法为：用二阶