CFD2015-第5讲-差分方法3

计算流体力学讲义2015第五讲差分方法（3）李新亮lixl；力学所主楼219；82543801,知识点：激波捕捉格式TVD、WENO、MUSCL、NND,1,CopyrightbyLiXinliang,课件下载：,知识回顾,1.差分格式的分辨率,有效网格点数：一个波长里面的网格点数（PPW:PointperWavelength),修正波数,2.群速度控制(GVC),CopyrightbyLiXinliang,3,3.流通矢量分裂,“在二阶迎风与二阶中心格式中选一个”“选接近一阶迎风的”,GVC2格式,CopyrightbyLiXinliang,4,5.1非物理振荡及TVD格式,1.数值解中的非物理振荡,间断附近非物理振荡的根源,理论1：色散误差导致各波传播速度不同（第4讲）,理论2：粘性耗散不足,思路：物理问题有粘；物理粘性足以克服本身振荡数值方法错误计算了物理粘性不足以克服振荡,物理问题本身也可能振荡。但如果错误计算物理粘性，则会错误地加剧（或衰减）振荡。,1）非物理振荡的原因分析,理论3：格式不能保单调,CopyrightbyLiXinliang,5,数值实验,二阶中心差分,计算域0,1,网格点201（Dx=0.005）时间步长Dx=0.0005,T=0.1时刻的u分布,Re=200Dx=0.005,现象：Dx一定时，减小Reynolds数可抑制振荡Reynolds数一定时，减小Dx可抑制振荡,暗示,是某一特征量,Re=2000Dx=0.005,Re=2000Dx=0.0005,相同,CopyrightbyLiXinliang,6,对流-扩散方程的特性：,n,n+1,（线性）差分方程：,某点的值是上一时刻周围几个点上值的线性组合,物理上要求系数ak均非负,含义：某处浓度的增加对下一时刻周围浓度的影响为正。,j-2j-1jj+1j+2,差分方程单调性（无振荡）条件：差分方程（1）中的系数非负,网格Reynolds数,CopyrightbyLiXinliang,7,2）重要概念：网格Reynolds数以网格尺度度量的Reynolds数,含义：数值振荡流动尺度为网格尺度网格Reynolds数小，该尺度的能量被耗散掉不发生振荡,j,j+1,j-1,过于苛刻的条件,单方向网格点数106，三维1018,单纯靠物理粘性抑制振荡，网格间距必须足够小，通常难以实现,网格足够小：不会发生振荡；网格小于激波的实际厚度，则不会振荡,网格Reynolds数足够小时，物理粘性发挥作用，抑制振荡,CopyrightbyLiXinliang,8,3）人工粘性,物理粘性足够小才发挥作用，Reynolds数很高时很难做到,思路：人为增加粘性系数（添加人工粘性）抑制振荡,优点：方法简便，有抑制振荡效果缺点：改变了物理问题，带来误差,湍流、分离流等对粘性敏感：非物理解,分离流对粘性敏感,转捩对粘性敏感,很难计算对粘性敏感的问题,改进措施：A:局部施加人工粘性B:高阶人工粘性,VonNeumann,MacCormack,人工粘性系数,大梯度区，加大人工粘性,光滑区为二阶小量,人工粘性项,CopyrightbyLiXinliang,9,4)数值振荡的定量描述总变差,对于离散函数uj定义总变差:,单调函数,振荡函数,j=1,j=N,含义：反映了振荡的剧烈程度,双曲型守恒方程,特点：沿特征线，u不变,特征线未相交总变差不变,特征线相交总变差减小,结论：单个双曲型方程，总变差不增（TotalVariationDiminishing：TVD）,CopyrightbyLiXinliang,10,2概念：单调格式、保单调格式与TVD格式,n时刻：单调函数,j=1,j=N,n+1时刻：仍是单调函数,j=1,j=N,设n时刻是单调的，如果n+1时刻的解仍保证单调，则称该格式为保单调格式。,保单调格式,基本结论：常系数的单调格式只能是一阶单调格式必是保单调的；线性格式，单调与保单调等价,格式：如果满足则称其为单调格式。,单调格式：,单调格式,保单调格式：,TVD格式,总变差不增,TVD,保单调,单调,CopyrightbyLiXinliang,11,3.TVD格式的理论基础Harten定理,Harten定理：,如果差分格式可写成如下形式：,且,则格式（1）是TVD格式,（1）,可验证：Roe格式是TVD格式,保证“系数非负”,含义：“单调格式必是TVD格式”,CopyrightbyLiXinliang,12,例7.2.1：,考虑线性单波方程：,试讨论如下Lax-Wendroff格式,二阶中心,人工粘性,是否满足Harten条件,常系数的单调格式只有一阶精度,对比条件：,不满足Harten条件,CopyrightbyLiXinliang,13,知识回顾：Lax-Wendroff格式,Taylor展开，写出修正方程,时-空二阶精度,巧妙添加人工粘性，不但克服了不稳定性，而且抵消了时间误差，提高了时间精度,类似方法：Beam-Warming格式,人工粘性,二阶精度迎风差分,人工粘性，且提高时间精度,特点：全离散、时刻耦合,CopyrightbyLiXinliang,14,4.构建TVD格式,思路：对现有格式进行改造，使之符合Harten条件,通常在Roe、L-W、B-M（或其组合）基础上改进80年代初、这些格式是主流,原格式（2阶）=1阶迎风+修正项,新格式=1阶迎风+限制函数*修正项,1.以二阶中心及二阶迎风格式为基础的改造,2阶迎风,2阶中心,2个候选格式：,思路1：两个里面选一个(GVC2,NND2)思路2：利用二者的组合,改造、新格式,二阶迎风,二阶中心,CopyrightbyLiXinliang,15,显然格式为2阶中心可验证：格式为2阶迎风,新格式：,二者组合仍为二阶,根据Harten定理，可知,时，可满足TVD性质,(2)精度条件,二阶精度区,TVD区,二阶精度TVD区（二者交集）,且,CopyrightbyLiXinliang,16,限制器（limiter),二阶中心,二阶迎风,GVC2格式,NND2格式,GVC2与NND2格式的限制器,CopyrightbyLiXinliang,17,常用的限制器,Minmod型；,Superbee型,VanLeer型,其他：如VanAlbada,CopyrightbyLiXinliang,18,概念：保单调区,1阶迎风的预测值,2阶中心的修正量,2阶迎风的修正量,精度差，但鲁棒性好,精度高，但有些情况下预测结果“不靠谱”,作为“标杆”检验高阶修正量是否可用,趋势相反时，不可用；相差超过2倍时，不可用,CopyrightbyLiXinliang,19,历史上，TVD格式是在Roe、L-W、B-M（或其组合）基础上改进80年代初、这些格式是主流,2以L-W格式为基础改造的格式,L-W,引入限制函数（限制器）,1阶迎风部分,修正项,CopyrightbyLiXinliang,20,显然格式为LW（2阶）可验证：格式为B-M（2阶）,新格式=1阶迎风+j*（LW格式-1阶迎风）,新格式：,LW,BM均为线性格式，二者组合仍为二阶,根据Harten定理，可知,时，可满足TVD性质,(2)精度条件,Beam-Warming,二阶精度区,TVD区,二阶精度TVD区（二者交集）,CopyrightbyLiXinliang,21,5.2其他激波捕捉格式简介,1.Godnov格式,思路：利用精确Riemann解难点：精确Riemann解要求间断两侧物理量为常数对策：采用分片常数代替原先的函数,方法：1）采用分片常数代替原先的函数即假设xj-1/2,xj+1/2区间物理量为xj点上的值2）在每个间断处，精确求解Riemann问题，得到Dt时刻物理量的分布。（假设Dt足够小，各间断之间无相互影响，因此可独立求解）。3）在区域xj-1/2,xj+1/2上平均，得到Dt时刻xj点上物理量的值。4）重复1-3，直到指定时刻。,j-1jj+1,间断1间断2间断3,2.NND格式,“在二阶迎风与二阶中心格式中选一个”“选接近一阶迎风的”如果二阶中心与二阶迎风给出的修正趋势相反，干脆不修正。,Minmod(a,b):a,b符号相同时，取绝对值小的；符号相反时，取0.,k=1/3三阶迎风,k=0,k=-1,CopyrightbyLiXinliang,23,2.MUSCL格式（VanLeer）,二阶迎风,Formm格式,k=1二阶中心,1阶迎风,修正部分,间断处都会产生振荡,思路：振荡由修正项引起，需要对修正项进行限制,CopyrightbyLiXinliang,24,新格式,函数minmod(a,b)a,b符号相反，则返回0符号相同，则返回绝对值小的,1阶迎风,修正部分,j-1,j,j+1,j-2,j+1/2,思路：修正部分是由j-1,j,j+1三点计算而得，是与的插值,旧格式,二者的组合：j-1,j计算；j,j+1计算,如果两个修正方案趋势（符号）相反，干脆不修正如果趋势相同，则取（绝对值）最小者避免振荡,两个修正方案,选前者：二阶迎风选后者：二阶中心组合之：Formm、3阶迎风,CopyrightbyLiXinliang,25,含VanAlbada限制器的3阶MUSCL格式,推荐方法：,VanAlbada限制器,光滑区相差不大,间断区，相差很大,光滑区相差不大，s趋近于1,间断区相差很大，s趋近于0,光滑区，趋近3阶迎风,间断区，趋近1阶迎风,理论上并不满足TVD性质；但实际效果很好,a0）,负通量情况（a0）的方法计算,采用针对负通量（a0）的方法计算,可采用Steger-Warming,L-F,VanLeer等分裂，见第4讲,1)计算,，它们是的函数（推荐使用Roe平均计算）,效果更好但计算量较大,计算量小但效果略差有轻微振荡,CopyrightbyLiXinliang,38,5.WENO格式的改进,MartinP的WENO-SYMBO格式（MartinP,JCP220,270-289,2006),思路：原WENO格式采用迎风型网格基，耗散偏大新的WENO格式采用对称型网格基，耗散小采用优化方法，提高格式的分辨率,理想权重情况下差分格式的修正波数；(SYMOO:未优化8阶中心差分；SYMBO分辨率优化后的4阶格式）,针对如下Sod激波管问题,用5阶WENO格式计算其数值解，画出t=0.14时刻密度、速度及压力的分布；并与精确解进行比较（要求数值解与精确解画在同一张图上，便于比较）。要求：空间网格数100，时间推进格式选用3阶Runge-Kutta,时间步长自选。,39,CopyrightbyLiXinliang,作业5.1,针对如下Shu-Osher激波-密度扰动波干扰问题：,用5阶WENO格式计算其数值解，画出t=0.1时刻密度、速度及压力的分布要求：1）空间网格数200，时间推进格式选用3阶Runge-Kutta,时间步长自选。2）使用2000个网格点计算，其结果作为“精确解”，与其它结果画在一起，便于比较。,40,CopyrightbyLiXinliang,作业5.2,初值为:,41,CopyrightbyLiXinliang,作业5.3（选作题）,推导7阶精度的WENO格式，给出详细的推导过程及格式的具体表达式,提示：与5阶WENO推导思路相同，但网格基架点扩大了，模板数目也增加了,j-3j-2j-1jj+1j+2j+3,正通量(