计算流体力学(中科院力学所)-第10讲-有限体积法2

计算流体力学讲义第十讲有限体积法（2）李新亮[e-mail]；力学所主楼219；82543801知识点：1讲义、课件上传至(流体中文网）-“流体论坛”-“CFD基础理论”讲课录像及讲义上传至网盘近似Riemann求解器——HLL/HLLC中心型有限体积法——人工粘性具体问题的计算——翼型绕流(边界条件、具体解法、粘性项处理)CopyrightbyLiXinliang2知识回顾0y(U)fx(U)ftU210F(U)tU在以某节点为中心的控制体上积分i,jk非结构网格的控制体i+1,ji-1,ji,j+1i,j-1k3k1k2k4k5结构网格的控制体0][dF(U)tU01dsnFtUdUU1dsdnFF(U)xyn01mmijHtU体积平均))()((1ymxmmmmnUnUss2ffnFH控制体边界垂直于节点连线（也可选其他方式）垂直平分线n1)建立控制体mx2)在控制体上积分离散方程重构：由节点上平均值给出函数分布，最终给出通量ijU表示第m个界面上的值1m2m3m4m1.有限体积法的离散过程——重构1）重构——两种不同的重构方案，向左偏及向右偏。给出两种结果：及CopyrightbyLiXinliang3i,ji+1,ji-1,ji,j+1i,j-1n左重构右重构LmU2）由左右重构得到的自变量：和给出通量方案A:FVS方案B:解Riemann问题（常用）LjiU,2/1RjiU,2/1LjU2/1RjU2/1例如：0阶重构：jiLjiUU,,2/1jiRjiUU,1,2/1线性重构:)(21,1,,,2/1jijijiLjiUUUU)(21,1,2,1,2/1jijijiRjiUUUU用i,i-1点的值插i+1/2点的值（网格剧烈变化时，应当用实际坐标插值）)()(2)(,1,,1,,1,,2/1jijijijijjijiLjiUUxxxxUU用i,i+1点的值插i+1/2点的值1111,,,pvu2222,,,pvuxy看似二维Riemann问题，其实是一维的，坐标旋转一下就行了RmU)))(ˆRmLmmU(Uf(Uff2.迎风型有限体积法（称为数值流通量）的含义CopyrightbyLiXinliang4重要概念澄清：重构与插值0)(xuftuA.有限差分法：xffxfjjj2/12/1ˆˆj+1/2切线j-1/2jj-12/1ˆjf)(2/1jxf2/1jf2/1ˆjf注意：与f在xj+1/2点的值含义不同！2/1ˆjf用周围几个点的值计算的过程称为“重构”，不能理解为用来插值2/1ˆjfjfjf)(2/1jxf记号确实容易混淆，让人容易联想起。记为更好些2/1jf)(2/1jxf2/1ˆjf否则，最高只能达到2阶精度了！是控制体内的平均值（称为数值流通量）的含义CopyrightbyLiXinliang5重要概念澄清：重构与插值0)(xuftuB.有限体积法：02/12/1xfftujjjj+1/2j-1/2)(2/12/1jjxff2/1jf2/1))((2/12/1jxxjjfxuff确实为f在xj+1/2点的值！通常做法：1）用计算出2）ju2/1ju)(2/12/1jjuff)(2/12/1jjxuuu在xj+1/2点的值！关键：是用计算（称为重构），而不是用计算（是标准的插值）；否则最高也只能达到2阶精度。ju2/1juju2/12/1)(1jjxxjdxxuxujuju1ju1ju2/1ju6概念：MUSCL与非MUSC类方法0)(xuftuxffxujjj2/12/1ˆˆj+1/2切线j-1/2j-12/1ˆjf2/1jf2/1ˆjfxffxujjj2/12/1差分有限体积juju方法1（非MUSCL类）：直接利用周围几个点的函数值或）直接计算（或）如何计算或？2/1ˆjf2/1jf方法2（MUSCL类）：利用周围几个点的自变量值（或）计算出（或）；然后再计算（或）jf)(juf2/1ˆjf2/1jfjuju2/1ˆju2/1ju)ˆ(ˆ2/12/1jjuff)(2/12/1jjuff当f=f(u)是连函数时，二者精度相同uCufuuff)()(00f的误差与u的误差同阶CopyrightbyLiXinliang710.1近似Riemann解迎风型有限体积法通常需要求解Riemann问题；精确Riemann解计算量大10.1.1HLLRiemann近似求解器（Harten,Lax&vanLeer）00)0,(xxxRLUUU0xtf(U)U思路：在控制体上积分，分析积分方程txRxLxRTZLTZT在控制体上积分Euler方程TxxRL,0,0)(0dxdtxtTxxRLf(U)U0))],()),([)]0,(),([0TRLxxdttxftxdxxUTxURL(Uf(URef.：E.F.Toro:RiemannSolversandNumericalMethodsforFluidDynamics,Springer,2009(ThirdEdition)控制体内质量（动量、能量）的减少等于流出控制面的通量LxRxCopyrightbyLiXinliang8txRxLxRTZLTZT0))],()),([)]0,(),([0TRLxxdttxftxdxxUTxURL(Uf(U若控制体空间足够大（或时间跨度足够小），扰动波未达到控制体的边界（如图）未扰动的值00)0,(xxxRLUUU边界上，未受扰动，保持常数0),(RLffUUUTTxxTxLLRRxxRL))();((RRLUffUffL把积分域分成三段；中间为受扰动段为扰动波传播的速度（激波或稀疏波的波头传播速度）RLZZ,0),(),(),(RLffUUUUUTTxxdxTxdxTxdxTxLLRRxTZTZTZTZxRRRLLL左、右波的速度假设已知LTZRTZLxRxCopyrightbyLiXinliang9txRxLxRTZLTZTtxRxLxRTZLTZT未扰动段，不影响积分，可不考虑0),(),(),(RLffUUUUUTTxxdxTxdxTxdxTxLLRRxTZTZTZTZxRRRLLL0)(),()(RLffUUUUUTTxxTZxdxTxxTZLLRRRRRTZTZLLLRL未扰动区，保持常值ULURLRLLRRTZTZLRZZTZZdxTxTZTZRL)(),()(1UUffURL消项TZZdxTxLLRRTZTZRL)(),(UUffURLT时刻，扰动区内的平均值T时刻的流动以平均值代替分布值CopyrightbyLiXinliang10思路：以扰动区内的平均值，代替瞬时值（扰动区很小，误差也不大）HLLU含义：T时刻扰动区内的平均值LRLLRRTZTZLRhllZZTZZdxTxTZTZRL)(),()(1UUffUURLRRRLhllLLZtxifZtxZifZtxiftx///),(~UUUU于是，HLL的近似Riemann解为：HLL近似解是左右波都是激波，中心区域均匀分布的一个近似模型txRxLxRTZLTZTtRxLxRTZLTZ双激波近似HLLULURU假设：1）左右波都是激波2）中心区无接触间断（近似）显然，假设中心区无接触间断是不合理的（造成方程数多于未知数个数），但可作为一种近似（类似最小二乘解）CopyrightbyLiXinliang11利用HLL近似Riemann解计算中心线上的通量xRxLxRTZLTZTdTtUT),0(00f在控制体（图中红色部分）积分Euler方程，得：TxL,00,0)(00dxdtxtTxLf(U)U0))],0()),([)]0,(),([00TLxdttftxdxxUTxUL(Uf(U00LLLLLhllTTTZTZffUULLhllLLZfUUf)(0控制体积（图中绿色部分）积分TxR,0,0RLhllRRZfUUf)(0二者计算结果相同带入表达式HLLULRLRLRLRhllZZZZZZ)(UUfffRLCopyrightbyLiXinliang12RLRRLRLRLRLRLhllZZifZifZZZZZZZifi00)(02/1fUUffffRLLLZRZ最终界面的HLL通量如下：原理：1）双激波近似，假设初始间断演化成两道激波，之间物理量为常数；2）利用积分关系式，确定了两道激波（左行激波与右行激波）之间的物理量；3）利用积分关系式，计算了穿过x=0界面的流通量优点：计算简单缺点：没有考虑到接触间断。该现象（初始间断产生两道激波、中间没有接触间断）不可能发生。Riemann问题是三波问题，该近似为两波问题，假设过强。CopyrightbyLiXinliang1310.1.2HLLCRiemann近似求解器（Toro）发展了HLL近似解，用三波模型来近似（如图）RZLZ*Z三波近似，左、中、右波的波速为，，中间波为接触间断LZRZ*ZLURUL*UR*URRRRLLLLZtxifZtxZifZtxZifZtxiftx////),(~****UUUUUHLLC近似解RRRRLLLLhllciZifZZifZZifZif0000****2/1fffffT时刻的流动状态LUL*UR*URUCopyrightbyLiXinliang14xRTZLTZ*ZRZLZ根据积分关系，可知L*UR*U)(**LLLLLZUUff红色区域积分可得)(**RRRRRZUUff蓝色区域积分可得假设中间波为接触间断：******uuupppRLRL**uZ利用积分关系计算接触间断的速度及其左右的物理量如果把ZL,ZR也当成未知数，可精确求解（6方程6未知数，Riemann精确解），但计算复杂。R-H关系式；弱解定义式含义：控制体内质量的增加等于控制体1控制体2CopyrightbyLiXinliang15问题最终描述：求解超定方程——4个未知数，6个方程RZLZ*ZLURUL*UR*U)(**LLLLLZUUff)(**RRRRRZUUff求解方程组：RLKKK,)**f(Uf其中：假设波速ZR,ZL已知，未知数——4个)(),(,,*********ZuuupppRLRLRL与精确Riemann解的区别：本近似解假设左、右波速度已知，因此，未知数只有4个；精确Riemann解认为ZL,ZR也是未知数，因此6个方程，6个未知数6个方程，4个未知数，怎么解？6个方程求解思路——只使用其中的4个方程方案（1）：利用的连续性方程及动量方程（共4个），解出全部4个未知数求解方案不能太复杂，复杂度不能超过求解原6个未知数的方程（Riemann精确解）；原则——务必给出解析解CopyrightbyLiXinliang16)(**LLLLLZUUff)(**RRRRRZUUff