水动力学计算薛雷平

水动动力学计学计算若干问题问题2012年研研究生暑期学学校上海交通大学学 薛雷平提纲数值计算要点有限体积法压力修正法大型稀疏方程组迭代求解12/9/202221 数值计算要点基本模型的选择计算方法的选择边界条件和计算流程的选择网格生成及其收敛性研究结果分析和总结12/9/202231 数值计算要点基本模型的选择基本模型的选择 可压缩非定场粘性流动方程不可压缩NS方程+湍流模型（RANS）势流方程（理想流体）边界层方程2D简化、几何简化、定解条件简化适用范围适用范围?简化模型的误差估计简化模型的误差估计?可通过典型算例不同层次的简化计算 得到合适的基本计算模型12/9/202241 数值计算要点计算方法的选择计算方法的选择离散方法：离散方法：有限体积有限体积，有限差分，有限元，有限差分，有限元计算格式：计算格式：伪可压缩模型的时间推进法（显式），伪可压缩模型的时间推进法（显式），压力修正法压力修正法（隐式）（隐式）湍流模型、多相流及空化模型、数值造波湍流模型、多相流及空化模型、数值造波消波消波12/9/202251 数值计算要点边界条件和计算流程的选择 边界位置及给出的边界条件并不一定要与实际问题完全一致 通过物理意义的考虑和试算得到能够简单实现，但又不使实际问题失真的边界条件 注意点：来流速度分布和湍流度 非远场出口边界条件 对称几何的非对称流动 空化流动时的压力条件的实现 计算流程：对于复杂问题，通常可以先计算物理、几何简化的初始解，逐步加入复杂性12/9/202261 数值计算要点网格生成 截断误差取决与于 网格间距 和 解的梯度，合理分布网格，通过初步计算重新分布网格通常是必要的 粘性流动计算时边界层网格（最好是结构化网格），不同湍流模型对y+的要求网格收敛性研究 网格收敛性研究指采用不断加密的网格使解收敛 不是指采用某一密度的网格时解与实验值一致 这是获得可靠数值计算结果必不可少的步骤12/9/202271 数值计算要点结果分析和总结 数值计算结果分析与实验结果分析类似 尽量采用无量刚表达式 CFD不是colour fluid dynamics，少用云图，多用曲线 分析总结的层次：分析总结的层次：1：初步的定性规律（相关性和大小）2：进一步的定性规律（线性、幂次、指数，先用量纲分析和理论储备获得基本规律，多变量时考虑采用变量分离法）3：定量规律（系数的拟合，加入特征点）4：分析定量规律背后的物理机理12/9/2022812/9/202292 有限体积法流体力学基本方程的一般形式 对流扩散方程随体变化随体变化=时间变化时间变化+对流对流=扩散项扩散项+体积源项体积源项+面积源余项面积源余项AVssutDtD12/9/2022102 有限体积法对流扩散方程的体积分形式 V为任意空间体积（控制体），F为体积V的表面该输运方程该输运方程可由 -微分方程体积分后直接由GAUSS积分公式得到，-也是控制体中物理量变化规律控制体中物理量变化规律的描述dAsndVsdAndAuundVtdVDtDAAVVAAAtVV)(12/9/2022112 有限体积法将积分方程应用于一个微元有限体积PV1a2a,vPPsffAPVffffAPsaVsauuatV12/9/2022122 有限体积法对流项的离散 迎风格式 +高阶修正faLPRP)0,min(),0,max(:ffffffAffRPfLPffAFmFmFmFmauuaFmFmFmFmuua的质量流量通过12/9/2022132 有限体积法扩散项的离散梯度格式faLPxRPx上位于逆变基协变基，局部坐标，fLPRPiiiffffffgxxggggggaaA,:3213322111g2g3g12/9/2022142 有限体积法扩散项的离散梯度格式agaxxFFFFFagagagaggggggggAgAggggggfRPLPRPLPRPfffffLPRPfffffffff1111111111113322111321321,-1,线性插值获得由12/9/2022152 有限体积法扩散项的离散梯度格式1111,gagaagaagaafffffffffLPRPfff时，12/9/2022162 有限体积法有限体梯度项的数值计算PfffPAVVadAndV12/9/2022172 有限体积法时间变化项的离散 )(243)()(2111211121tOtttOtttOttnPnPnPnPnPnPnPnPnPnPtVtVtVPPPPP2 有限体积法时间变化项的离散与对流项合并（对任意f:LP=P)12/9/202218fffPffPPPFmFmFmtV-fffffPRPfPPffAPFmFmtVuuatV归纳可用于结构和非结构网格2 有限体积法12/9/202219迭代中间解的线性方程组，构成一阶时间向后）111111111-(nnffnffffffnfffnPnRPfffnPnPnPffffAPgaFmFmtVauuatVfPfnPPPfffPfPfnfRPPfnPPPAtVAFmArAA111,-，）（PV1a2aPLP,)1(fRP)2(fRP3 压力修正法SIMPLE算法全名为压力耦合方程组的半隐式方法（Semi-Implicit Method for Pressure Linked Equations）被广泛使用的求解不可压缩流场的数值方法，于1972年由S.Patankar与B.Spalding提出原因：流体不可压缩假设使得扰动传播速度趋于无穷大，方程数学性质由双曲型变为椭圆型，根据特征线方法设计的气动计算方法无法直接使用特点：必须求解压力修正量的隐式方程，每个时间步必须多次迭代计算动量方程、压力修正方程及其他辅助方程；可直接进行定常计算（时间步长无穷大）12/9/2022203 压力修正法来由1 速度由NS方程求得2 压力必须由连续性方程求得代入NS方程的p,必须使得解出的u满足连续性方程12/9/202221howpuuupuutu?013 压力修正法基本步骤0 假设初始解p*1 代入NS方程求得速度初始解u*2 获得速度解与压力解的离散关系3 真解为初始解与修正解的叠加4 根据2：设12/9/202222PPffRPPfPPPrpVuAuA）（-uuuppp,SIMPLECSIMPLEAAApVuAfPfPPPPPPPPP:1,:0,5 该假设同样适用于控制面 f6 代入连续性方程若 ，跳出迭代ffpfPfffPffppAVupVuAffffpfPfffffffffffuaAVpauauauua*0*0*fffua3 压力修正法7 回顾 舍去非正交项可得 压力修正方程为12/9/2022231111,gagaagaagaafffffffffLPRPfff时，PRPffpfPfpfPffppAVAVpafPfPPfpfPfPfffffRPPfPPPBBAVBuapBpB,*-）（3 压力修正法8 隐式求得p，修正压力9 修正速度10 将p作为新的初始解p*，回步骤1备注：使用非交错网格时，计算控制面f的速度时，要考虑压力梯度的作用，否则会产生振荡解12/9/202224refPPPPppppPPPPPPpAVuu*fffPfPffffffPfPfPRPRPRPPRPRPPPPPPppAVuuHHpAVuHpAVuHpAVu,可得：3 压力修正法控制面体积流量计算归纳为12/9/202225fPRPffPfPffffffPfPfffffffffPRPfffpppgAVuappaAVuauaFpgpapppa114大型稀疏方程组迭代求解线性方程组的迭代求解：采用迭代矩阵进行求解：12/9/202226rxAxxxxArxMkrxMAxMrxMAMMMAMAkkkkk11-,，或为迭代次数为迭代矩阵4大型稀疏方程组迭代求解对M的要求：容易求逆M与A尽量接近 或 M*x与A*x 尽量接近Jacobi：M为A的主对角线Gauss-Seidl：M为A的主对角线+下三角元素12/9/202227xAxMAM或4大型稀疏方程组迭代求解非完全三角分裂法和SIPSIP:strongly implicit procedure，Stone1968针对结构网格上形成的5对角线（2D）或7对角线（2D）系数矩阵的高效迭代方法非完全上下三角分裂法的拓展12/9/2022284大型稀疏方程组迭代求解5对角线系数矩阵4大型稀疏方程组迭代求解非完全上下三角矩阵分裂 M=L U=A+N ;0;0,0,0ijijijijijijAifAMAifULSIP:要求 M A (A)P=且设 对比系数有 4大型稀疏方程组迭代求解WNWWESEENNWNSSESPSENWPAMMAMMAMMAMMAMMMPANW N W P E S SE4大型稀疏方程组迭代求解代入4大型稀疏方程组迭代求解可得=0 时，为普通的非完全三角分裂法4大型稀疏方程组迭代求解SIP方法的迭代矩阵M，可结合共轭梯度类方法结合，作为CG法（对称系数矩阵）以及BiCGSTAB（非对称系数矩阵）的前置优化矩阵使用。尤其对于多块结构化网格，可以取得更快的迭代收敛速度。12/9/202235 谢谢！12/9/202236BICGSTAB主要思想12/9/20223712/9/202238BICGSTAB主要思想BICGSTAB计算流程12/9/202239BICGSTAB计算流程12/9/202240CG 计算流程12/9/202241