水利机械计算水动力学------提交报告

1、一、求解N-S方程的数值方法求解 N-S 方程的数值方法可分为有网格和无网格两类。有网格方法的求解思路是将连续空间划分为离散网格 (或单元)，并在网格点 (或单元) 上离散控制方程，从而寻求控制方程的数值解；无网格方法则将着眼点回归到流体运动的基本单元质点上，通过追踪质点运动性质的变化过程来获得流体解答。网格数值模拟方法发展较早，求解方法较成熟，应用也较广；而无网格方法则刚刚起步，具有很大的发展空间。下面分别就两类中较有代表性的方法予以介绍。1.1 网格数值模拟方法数值计算网格可分为不随流体运动而改变形状的欧拉网格和随流体运动改变形状的拉格朗日网格两类。欧拉网格划分较简单，但在用于具有自由表面

2、流体运动数值计算中需要对自由表面位置进行适时追踪；拉格朗日网格虽然无需进行自由表面追踪，但需要在计算过程中随时调整计算网格，以适应流体自由表面形状，额外增加了计算量，并且在流体变形剧烈处还会出现网格失效现象。有限差分法 (Finite Difference Method, FDM)和边界元法 (Boundary Element Method, BEM) 都是欧拉网格方法中常用的方法。在拉格朗日网格法中，有限元法和边界元法可直接应用，而有限差分法一般需对坐标做适体变换。用欧拉网格法求解含有自由表面的水波动力问题难点在于追踪自由表面。目前比较有代表性的自由表面追踪方法有 PIC (Particle

3、-In-Cell) 法、MAC (Marker-And-Cell) 法和VOF (VolumeOfFluid) 法等。以下分别予以介绍。1.1.1 质点网格法 (Particle-in-Cell，简称 PIC 法) PIC法是 Harlow提出的一种用欧拉矩形网格计算多种介质流体运动的方法。该方法把流体既视为连续介质，又视为带有一定质量的质点，然后研究质点在经过固定欧拉网格上的运动性质。PIC 法具有计算多相流和处理三维自由表面的能力，曾经成功运用于模拟二维流体发生剧烈变形的情况。其缺点是占用较多计算机内存；模式只能给出自由表面单元位置，而不能给出自由表面精确位置；另外在流体变形剧烈处易出现较

4、大插值误差。1.1.2 MAC 法 (Marker-And-Cell) 针对 PIC 法的上述缺点，20 世纪 60 年代中期，Harlow, Welch提出了一种改进方法，称为 MAC 法。最初的 MAC 法在所有单元内部都布满无质量的标记点，通过跟踪这些点即可判断自由表面的位置。MAC 法在求解时以速度分量和压力为因变量，在固定网格上采用有限差分法离散控制方程。在 MAC 法之后，人们又先后提出了 SMAC 法、ABAMC 法、SUMAC 法和 TUMMAC 法11系列。新版本模型改为只在自由表面处设标记点，对这些点的连续跟踪可得到自由表面的准确位置，且所耗机时大为减少。Hirt 和 Ni

5、chlos在总结当时的自由表面跟踪方法的基础上还提出了“线段法”(Line Segment Method) 的概念，可近似确定自由表面的位置，但这种方法推广到三维的情况较困难。MAC 法于 20世纪 60 年代提出后先是用于涌潮的传播研究，后自 80 年代始大量用于波浪研究。如，Miyata运用一种二维 MAC 格式研究了波浪破碎；Sakai等也提出一种基于早期 MAC 法的数值模型，来模拟水波泼溅以及波浪破碎过程中二阶乃至三阶漩涡的生成。Gao 和 Zhao运用二维 MAC 模型，研究了波浪、建筑物及沙滩的交互作用。在 MAC 法之后，学者们开发了更具通用性和计算效率的模型。如，Chorin

6、提出的投影法 (Projection Method) 就是以其数学上的精确、简明而著称的一个。以 MAC 法为基础，Hirt 等和 Nichols 等开发了易于操作的 SOLA (solution-algorithm) 模型。Kothe等利用不完全 Cholesky共轭梯度法，将表面张力模拟为一种体积力，开发出了 RIPPLE模型。1.1.3 VOF 法（Volume-Of-Fluid）1.1.3.1 算法思想Nichols 等是最早用利用VOF法的研究群体。VOF 法假设在整个流动区域内流体密度为一常数 ，空气密度设为 0，引进流体体积参数 (或者是 VOF 函数) ,从连续方程出发建立关于

7、F的对流方程：(1-1)F 在流体单元内为 1，在空单元中则为 0，而在自由表面单元中， F 值介于 0和 1 之间。给定速度场就可以通过上述对流方程确定F 值，进而确定任一时刻的自由表面单元位置。1.1.3.2 方法介绍 供体受体法 (Donor-Acceptor Method)：由于 VOF 函数在自由表面附近有很大梯度，用一般的激波捕捉方法虽可求解上述 VOF 对流方程，但由于存在严重数值耗散，难于精确确定自由表面位置。为此，Hirt 和 Nichols提出了供体受体法 (Donor-Acceptor Method)。该法首先依据流速在 x , y 方向的梯度大小确定自由水面的走向，并根

8、据流速的方向判断自由水面附近的单元格为供体格还是受体格，并进一步在运动方程求解中决定采用逆风或顺风格式。这就有效地控制了数值耗散。有关数值实现的详细信息可参阅 Hirt和 Nichlos、Kothe 、Liu 和 Lin等相关文献。倾斜界面法 (Sloping Interface Method)：上述供体受体方法对自由表面几何形状进行了简化处理，导致供体受体法精度不高。Youngs引入了倾斜界面法，较之于供体受体法，此类方法可更精确确定自由表面的位置。但目前在算法实现以及向三维扩展时有一定困难，也就妨碍了倾斜界面法的实际应用。层面设置法 (Level Set Method)：该法由 Sussm

9、an 等于 1994 年提出，是早期借助移动网格跟踪界面方法的改进版本。该方法很容易确定分界面并且数值耗散较小。但由于这种格式没有明确利用质量守恒方程，因此复杂情况下容易造成质量亏损。1.1.3.3 应用情况 Delft 水力学研究小组自 20世纪 90年代早期就利用 VOF 方法进行波浪动力学研究。Van der Meer 等提出一种模型，研究作用在沿岸建筑物上的卷破波。Iwata 等研究了由于水下建筑物所致的波浪破碎变形。Lin 和 Liu k - e 利用紊流模型封闭 N-S方程，开发了一种基于 VOF 法的模型，对倾斜海滩上孤立波破碎等诸多水波动问题进行了研究。Sussman等耦合 V

10、OF 法和层面设置法，用于三维情况的不可压缩流体的研究。上述方法虽然解决了自由表面的追踪问题，但当自由表面有剧烈变形时，各类方法在处理对流项时，数值耗散问题仍然都是不可避免的。1.2 无网格数值模拟方法无网格数值方法也称为粒子方法，该方法完全摒弃了数值网格的概念，由于运动方程中不直接出现对流项，自然也就不会出现网格方法中的数值耗散现象。目前文献上报道的粒子方法有 PAF法 (Particle-And-Force) 、无网格 Euler/N-S 方程解法、自由元 (Free Element) Galekin 法、离散元 (Discrete Element) 法、SPH(Smoothed- Par

11、ticle-Hydrodynamics) 法以及 MPS 法 (Moving Particle Semi-implicit Method)等。在水动力学研究中开展较多的方法主要有 MPS 和 SPH 两种，以下对其予以介绍。1.2.1 SPH 方法 SPH 方法的基础是插值理论。SPH 通过某个插值函数在质点上给出变量的近似值，从而将以偏微分方程形式给出的连续流体的守恒方程转换成积分方程。控制方程中的压力梯度项、扩散项等空间导数项也通过核子函数表示。Cummins 和 Rudman 提出了一种 SPH投影法 (SPH projection method)，通过求解压力泊松方程确定压力，从而确立

12、了具有严格算法的不可压缩模型。SPH 方法源于天体物理学，后扩展至水动力学研究。随后，Edmond 等用不可压缩 SPH模型求解 Lagrange形式的 N-S 方程，并将模型与大涡数值模拟相耦合，以模拟近岸孤立波的作用机理。继而，Shao 等将不可压缩 SPH 模型用于模拟带有自由表面的牛顿和非牛顿流体的研究。1.2.2 MPS方法 Koshizuka 等提出了一种粒子方法，用于计算不可压缩黏性流体，并将之命名为移动粒子半隐式方法 (Moving-Particle Semi-implicit method)，简称为 MPS 法。该方法中，粒子间相互作用通过一个核子函数来表征，梯度和扩散算子运

13、算都通过核子函数表达。流体的不可压缩性是通过控制粒子密度数不变实现的。计算过程分为显式 (Explicit Stage) 和半隐式两个阶段。显式阶段仅考虑重力和黏滞项，半隐式阶段用隐式方法求解压力泊松方程，但速度的二次校正仍然为显式。目前，该方法已用于波浪破碎机理和两相流的研究中。另外，Gotoh 还将亚格子尺度紊流模型引入 MPS 模型，获得了 N-S方程的 Lagrange 解。总体来看，水波动力学中的无网格数值方法还处于初期探索阶段，运用于实际工程中的例子还少见报道。但由于其独有的可很方便地追踪自由表面的特性，决定了其在某些特定的工程实际问题中有网格方法无可比拟的优越性，因而也具有很好的

14、发展前景。二、湍流及其数值模拟方法2.1 湍流湍流是流体的一种流动状态。当流速很小时，流体分层流动，互不混合，称为层流，也称为稳流或片流；逐渐增加流速，流体的流线开始出现波浪状的摆动，摆动的频率及振幅随流速的增加而增加，此种流况称为过渡流；当流速增加到很大时，流线不再清楚可辨，流场中有许多小漩涡，层流被破坏，相邻流层间不但有滑动，还有混合。这时的流体作不规则运动，有垂直于流管轴线方向的分速度产生，这种运动称为湍流，又称为乱流、扰流或紊流。2.2 湍流的数值模拟工程中绝大多数流体都是湍流。深层次地揭示流体运动的紊动结构及能量耗散过程是计算流体力学应用研究的极其重要的课题。出于不同的研究目的，湍流

15、数值模拟有以下三个不同的层次：直接数值模拟 ( Direct Numerical Simulation，简称 DNS)、雷诺平均数值模拟(Reynolds Averaged Navier-Stokes，简称 RANS) 和大涡数值模拟 (Large Eddy Simulation，简称 LES)。以下分别简要介绍。2.2.1 直接数值模拟由于 N-S 方程本身是封闭的，故从原则上讲可以求解所有湍流问题。用 DNS 法直接求解 N-S 方程，能获得最精细的流场信息。但在目前计算机发展水平下，进行 DNS 的应用研究是不现实的，目前只限于小规模的低雷诺数简单湍流物理机制的研究。DNS 常用的数值方

16、法是谱方法或伪谱法。Orszag 等最早用 DNS 计算了各向同性湍流。2.2.2 雷诺平均数值模拟由于雷诺平均方程是不封闭的，故 RANS 的核心思想是建立雷诺应力封闭模型，使得平均运动方程可解。目前比较常用的模型有零方程模型、一方程模型 (k方程模型)、二方程模型 (涡黏性模型，k -e 模型)、代数应力模型 (k -e - A模型)。其中k -e 模型是目前应用最广泛的湍流模型。在可以预见的将来，即使 DNS求解复杂湍流成为现实 (耗时很长)，能够快速地算出满足一定工程精度要求的湍流统计模式仍会受到工程师们的青睐。事实上，RANS 是目前工程界处理湍流问题的唯一方法。2.2.3 大涡数值

17、模拟LES 的基本思想是：把湍流瞬时运动通过某种滤波方法分解为大尺度运动和小尺度运动两部分。大涡运动通过直接求解 N-S 方程计算；小涡运动的影响概化为亚格子雷诺应力，需通过建立模型求解，模型称为亚格子尺度模型 (Subgrid Scale Model)。气象学家 Deardoff 首次把大涡模拟用于有工程意义的槽道中的流体运动数值计算中。1972 年起，Stanford 大学的 Ferziger和 Reynolds 领导的集体开始对 LES 做深入系统的研究。苏铭德 1982 年以来提出了一种代数应力模型，计算了槽道中的流体流动。Watanabe 用大涡模拟计算了三维波浪破碎。LES 至今在

18、气体动力学中的研究开展较多，在水动力学中的研究还十分有限，最主要的困难还不在于计算机的限制，而在其方法本身，如现有的亚格子尺度模型仍很不完善以及近壁模型的入流出流边界等问题。三、水力机械空化数值模拟问题描述：高速离心泵在工作时，容易在叶轮叶片处产生汽蚀，会严重影响离心泵的性能。在叶轮入口处配置诱导轮,可以提高离心泵的抗汽蚀性能。3.1 两相流基本理论在低比转数离心泵发生汽蚀时，从气泡的产生到消失时间极短，是一个复杂的动态气液两相流动过程。气液两相流流场的组成，从宏观看，包括气相、液相和气液界面3部分。气液界面的存在对两相流的流动特性有很大影响，流场中的界面在数量上和分布上经常发生变化。 在等温

19、两相流中，随流动的进行，压力呈现变化，两相界面的大小和分布状况即发生变化。当压力变化较小时，沿流道压降与介质压强相比，数值很小，界面及其分布可以看成是不变的。但在压力变化很大的两相流中，介质压力大小和气液两相界面的变化均很显著，因此流动特性也会随之改变。气液两相流的数值模拟包括气相场和气液间的相互干扰计算。相间干扰即气液两相间的动量、能量、质量的交换过程。常见算法可以分为欧拉一欧拉型算法和欧拉一拉格朗日算法两种类型。3.2 两相流控制方程连续方程：(3-1)式中：液滴蒸发进入气相的质量动量方程：(3-2)式中：F离散相粒子对气流的作用力3.3 气泡项体积分数方程(3-3)其中： 3.4 变螺距诱导轮的设计3.4.1 诱导轮提高泵抗汽蚀性能的原因高抗汽蚀性能的诱导轮在一定的程度上是以牺牲能量指标换来的，其效率并不高，但诱导轮本身负荷不大，可以改变叶轮进口流动情况，不明显降低水泵的效率。(a)无诱导轮 (b)有诱导轮图1 叶轮进口速度三角形图2 诱导轮模型由上图可知，加诱导轮之后，主叶轮叶片进口前相对速度和无诱导轮相比减小了。则主叶轮的汽蚀余量减小，泵的抗汽蚀性能提高。同时加诱导轮之后主叶轮叶片进口之前液体的绝对速度和未加诱导轮时相比也增加了，此的增加由诱导轮作用产生，并不引起压力下降。因此可以认为不影响主叶轮内最低