中科院计算流体力学最新讲义CFD114讲差分方法.ppt

《中科院计算流体力学最新讲义CFD114讲差分方法.ppt》由会员分享，可在线阅读，更多相关《中科院计算流体力学最新讲义CFD114讲差分方法.ppt（37页珍藏版）》请在三一办公上搜索。

1、计算流体力学讲义2011 第四讲 有限差分法（2）李新亮；力学所主楼219；82543801,知识点：离散误差的Fourier分析;间断周围数值振荡的原因、群速度控制（GVC）；模型方程向N-S方程的推广：流通矢量分裂,1,讲义、课件上传至(流体中文网）-“流体论坛”-“CFD基础理论”下载地址2：,Copyright by Li Xinliang,Copyright by Li Xinliang,2,知识回顾,1.有限差分基本原理,差商-微商,待定系数法,2.基本概念,差分格式、截断误差、精度,差分方程、修正方程,显格式、隐格式,全离散、半离散,守恒性、非守恒型,3.相容性、稳定性;Lax等

2、价定理,Copyright by Li Xinliang,3,4.1 差分格式的误差分析,概念澄清1：精度 vs.分辨率,描述 网格充分密集时差分格式的误差特性,描述网格有限的情况下，差分格式的误差特性,2 阶精度,误差,4阶精度,格式1,格式2,3阶精度,显然：足够小的情况下，格式1误差更小 并非足够小的情况下，格式2 有可能误差更小,精度特性,分辨率特性,Copyright by Li Xinliang,4,概念澄清2：有效网格,目的：计算差分，,Case 1,Case 2,网格点数增加了一倍，但问题也复杂了一倍,用同一差分格式计算，两种工况误差完全相同。,有效网格点数：一个波长里面的网格

3、点数（PPW:Point per Wavelength),有效波数,10个点,20个点,1.耗散与色散误差,5,Copyright by Li Xinliang,精确解,1阶迎风,2阶迎风,数值实验,时间推进：3步TVD型Runge-Kutta,且时间步长足够小（误差忽略）空间离散：1阶及2阶迎风格式（20个网格点）,实验观察到的现象 两类误差：振幅误差 相位误差（波速误差）,Copyright by Li Xinliang,6,对以上“实验现象”进行理论分析,半离散分析：假设时间推进是精确的，仅分析空间离散带来的误差（难度小、常用）全离散分析：同时分析时、空离散的误差（难度大）,考查问题：,

4、实际上就是普通三角函数，采用复数形式仅仅是为了理论推导方便。用实数形式 sin(kx),cos(kx)推导形式上略显繁琐。,精确解：,差分格式：,（1）,其他格式,假设对于：,有,隐含假设:线性差分格式 非线性系统作用于单波，会产生多个谐波,（2）,差分没有误差,Copyright by Li Xinliang,7,令：,（1）式化为：,“半离散化”：空间导数差分计算，时间方程（常微）精确计算,如果，无误差,分析（修正波数）与误差的关系,理想情况：的误差导致解的幅值误差 耗散误差 的误差导致解传播速度的误差 色散误差,假设对于：,有,反映了一个波内的点数。PPW（波内的点数）=,Copyrig

5、ht by Li Xinliang,8,耗散、色散误差分别由修正波数的实部和虚部决定。,关键参数：修正波数,含义：反应波数（谱）空间内差分的误差,任意函数：,定义：,求导数，精确解,差分解,Fourier 分析的任务计算出，并考差其与 的逼近程度。,考察格式分辨率（resolution）的重要指标 精度：反映 时的情况 分辨率：网格点数很少（例如波里面只有6个点）时的性能对于多尺度问题，分辨率更重要。牺牲精度，提高分辨率,优秀的差分格式，1个波长里面6个点 即可,精度 分辨率,Copyright by Li Xinliang,9,如何计算修正波数？,定义：,方法1.理论计算 根据差分具体表达式

6、及定义计算,例1：,令,则：,于是：,1阶迎风,例2：,2阶迎风,Copyright by Li Xinliang,10,方法2：数值计算,定义：,Step 1）选取计算域0,2p,计算网格（例如64，128）Step 2）给定波数 k,生成函数值Step 3)调用差分子程序，得到导数值Step 4)通过Fourier反变换，得到谱：,假设已有求差分的子程序（黑箱，已知是线性的）,线性黑箱,强调：研究CFD本身，不能只使用理论手段，还要用数值手段,根据修正波数的定义，有,Step 5)改变k的值，重复2-5，得到 对于 的依赖关系。画图,非线性情况会产生高次谐波，造成 step 4 中隐含的假

7、设无法成立 将Fourier分析手段拓展到非线性系统 需要研究的课题,隐含条件：只有波数为k的那个谱不为0（线性系统）,Copyright by Li Xinliang,11,中心差分格式的色散特性0：精确解；1：4阶普通2：6阶普通；3：4阶紧致4：6阶紧致；5：6阶超紧致,迎风差分格式的色散特性0：精确解，1：2阶迎风2：5阶迎风偏心 3：3阶迎风紧致4：5阶迎风紧致,每个波长里面2个网格点，谱方法的分辨率，差分法分辨率的极限（只有无穷阶精度才能达到）,20阶超紧致格式 接近谱方法,Copyright by Li Xinliang,12,不同差分格式的色散误差曲线,结论：要求分辨率相同的情

8、况下，采用高阶格式可放宽空间网格步长，从而减少计算量,重要方向：高分辨率差分格式,0：精确解1：2阶迎风2:3阶迎风3：3阶迎风紧致4：5阶迎风紧致,指定误差要求的情况下，不同差分格式能模拟的最大a(a 越大，所需网格越少）,附录：部分差分格式,表中的迎风差分格式均针对 a0,当a0时，需把下标的“j+k”换成“j-k”（例如把j+2 换成j-2,把j-1换成j+1);并在表达式前加上“-”号。例：,迎风偏斜格式：上游的基架点更多些（或上游权重更大）,13,Copyright by Li Xinliang,4.2 数值解的群速度及间断处数值振荡来源,对于：,有,修正波数,数值解,色散误差：数值

9、解传播的速度与精确解不一致,数值解传播偏快,数值解传播偏慢,0：精确解；1：2阶迎风；2：5阶迎风偏心 3：3阶迎风紧致；4：5阶迎风紧致,快格式(FST)：慢格式(SLW):混合格式(MXD):,特点：波数越高，误差越严重,1.色散误差与群速度,14,Copyright by Li Xinliang,t=0.5时刻的精确解及数值解,空间离散：五种不同格式；时间推进：3阶Runge-Kutta,【数值实验】波的传播问题,观察现象：1）高波数成分误差严重，低波数成分误差不明显；2）二阶Pade格式的解传播速度快于精确解，其余格式偏慢；3）迎风型格式有耗散，尤其是二阶迎风格式；,概念：群速度 波包

10、传播的速度,15,Copyright by Li Xinliang,2.间断附近数值振荡的来源,【数值实验】间断的传播,计算域0,1;计算网格点100 时间推进：3阶Runge-Kutta 空间离散：1）二阶中心差分,2阶迎风及2阶中心格式的色散特性,2）二阶迎风差分,16,Copyright by Li Xinliang,2阶迎风,精确解,2阶中心,过激波数值振荡的根源 色散误差导致群速度不一致,快格式,慢格式,波前振荡,波后振荡,=+,群速度控制的基本思路（群速度控制 GVC:Fu&Ma）：间断前、后分别采用快格式和慢格式，可有效抑制振荡 Zhuang&Zhang:抑制波动原则,示意图：间

11、断的Fourier分解,好思路,17,Copyright by Li Xinliang,利用GVC的思想构造可计算间断的差分方法,1）间断的前后判据,简易方法：则j点在间断左侧,j-1,j,j+1,则j点在间断右侧,原理：越靠近间断，振荡越剧烈,（a0 时，右侧为“前”）,2）根据GVC的思想构造格式,间断前：快格式；间断后：慢格式；,格式 GVC2,3）改写成为守恒型,非线性情况，通常守恒型效果更好,NND,格式 GVC2a,18,Copyright by Li Xinliang,4）a 0 时，同样思路构造（利用对称性，仅需把下标j+k换成j-k即可）,采用GVC2a(NND 2a)格式的

12、计算结果 消除振荡,使用NND2a（守恒形式）；NND2（普通形式）及1阶迎风格式的计算结果,将j换成j+1,19,Copyright by Li Xinliang,4.3 通量分裂技术,格式F+格式F-,对流项：信息（波）从上游传至下游 上游信息更重要 迎风差分扩散项：信息从中心向周围扩散 不区分上、下游 中心差分迎风差分优点：有效利用信息传播的方向，增强稳定性,微分与差分方程的影响域,N-S方程：,单波方程：,单波方程 一个波，容易判断波传播方向N-S对流项（Euler）方程组：多波问题，复杂,双曲方程组的原则 特征分解，找到独立传播的波,常系数矩阵A的情况 完全解耦，独立求解,变系数矩阵

13、A的情况 局部讨论,20,Copyright by Li Xinliang,Copyright by Li Xinliang,21,Step 1 针对模型方程构造差分格式,Step 2 将对流项分裂成正、负通量,方法1：流通矢量分裂(FVS),特点：对流通矢量 f(U)本身进行分裂,Steger-Warming,van Leer,Lax-Friedrichs,利用通量分裂技术，将模型方程推广到Euler(或N-S)方程,Copyright by Li Xinliang,22,Step 3:利用step1构造的格式,将正负通量离散化,方法2：通量差分分裂（FDS）,特点：对流通矢量 f 的导数进

14、行分裂,Godnov,Roe,HLL,HLLC,方法3：AUSM 方法,van Leer分裂法+压力项单独处理,1）矩阵分裂法,2）Riemann解法,利用差分表达式，计算,求解Riemann问题，获得通量,1.Jacobian 系数矩阵及其性质,重要性质,特点：A 可以像常数一样，和求导运算交换,23,Copyright by Li Xinliang,2流通矢量分裂(FVS),利用性质,=,+,优点：耗散小缺点：导数间断,方式 A:,特点：不必进行矩阵运算，计算量小,Steger-Warming 分裂,A:Steger-Warming 分裂,24,Copyright by Li Xinlia

15、ng,Steger-Warming 具体步骤（以一维为例）,已知,1）计算2）计算3）计算4）带入（1）式得到5）利用不同的迎风格式，分别计算,（1）,（后差，前差）,6）,计算,7）时间推进,25,Copyright by Li Xinliang,二维问题的steger-Warming 分裂,令：,则：,具体使用步骤，以计算 为例,令 计算特征值 分裂特征值，计算 带入左式，计算正、负流通矢量 计算,计算 设置，并注意,对于曲线坐标系,仅需令,三维问题同样处理 二维、三维具体 公式见傅德薰等计算空气动力学 4.7节（158-162）书中公式有一定的排版错误，使用前务必重新仔细推导！,26,C

16、opyright by Li Xinliang,B:Lax-Friedrichs(L-F)分裂,特点：,正特征值 负特征值,=,+,缺点：耗散偏大,局部L-F分裂，每个点上计算 全局L-F分裂，全局（一维）上计算,足够大,数学性质（光滑性）最好，但耗散偏大,常数,与迎风格式结合，等价于人工粘性,例如，可取,27,Copyright by Li Xinliang,方式很多,=,+,S-W:,L-F:,=,+,Van Leer:,=,+,28,Copyright by Li Xinliang,分裂后 失去了A的性质（可以像常数一样与求导交换）,FVS分裂：优点：无需矩阵运算，计算量小 缺点：分裂后

17、改变了特征方向，耗散大,利用了性质,一般情况下：,变系数，不能与导数交换,实质：没有做到解耦；只是把原变量重新组合，组合后波的传播方向的保证 f+向正向传播，f-向负向传播 缺点：由于未解耦，各变量的误差会相互传递,29,Copyright by Li Xinliang,概念澄清：流通矢量分裂本身不带来耗散，但其会影响到差分的耗散;,举例：,分裂过程,耗散,如果差分格式无耗散（例如都用中心差分），则通量分裂不带来耗散。,=,+,向上平移,向下平移,分裂,差分格式,耗散,分裂后的流场越偏离原先流场，则总体耗散越大,精确满足，不引入误差！,如使用低精差分度格式，则对分裂形式敏感（推荐使用特征分裂）

18、如使用高精度格式（低耗散），则对分裂形式不敏感（可使用逐点分裂）,30,Copyright by Li Xinliang,3.特征重构方法,常系数方程组：,完全解耦,变系数情况 局部冻结系数,在基架点上系数 不变,计算：,在差分基架点上Aj 不变，可按常矩阵处理,局部冻结系数,分别采用后差和前差,优点：严格保证（局部）特征方向，数值解质量好；缺点：大量矩阵运算，计算量大。,31,Copyright by Li Xinliang,通常写成守恒型差分，计算,在基架点上系数 不变,具体步骤：假设已知 U，且针对模型方程（线性单波方程）已构造出差分格式,（1）,1）计算出,教材130页的公式(6.1.

19、11-6.1.13),式中用到各变量在j+1/2的值（例如）可使用j,j+1 点值的算术平均（如）或Roe平均（教材6.4节）；由 计算；方法很多，例如前面介绍的 或,32,Copyright by Li Xinliang,均可,2）在网格基上计算,计算fj+1/2用到的点,注意，在该网格基上（例如k=j-1,j,j+1）保持不变,例如：,3)利用已构造好的差分格式，计算通量4）得到总通量,5）计算差分（j点处）,步骤的算法描述（注意：实际上是两重循环）do j=1,N do k=j-1,j+1(网格基，可以是更多或更少点)enddo enddo do j=1,N enddo,需要多次矩阵运算

20、，计算量大 守恒性好，耗散小，数值解质量好,33,Copyright by Li Xinliang,作业题1：构造高分辨率差分格式，并进行理论分析及数值实验,针对单波方程：,对于空间导数，构造出一种不超过6点格式；并进行Fourier误差分析，画出kr,ki的曲线。要求：精度不限；网格基架点数不超过6个；能够分辨的波数范围尽量宽；（即kr,ki曲线近可能接近准确解）给出差分的具体表达式，画出kr,ki的曲线；说明构造格式的阶数，并采用本PPT第5页的方法给出的精度验证；,形如：,另外，进行如下数值验证：,空间采用20个网格点，采用新构造的差分格式离散；时间推进采用3步Runge-Kutta方法

21、，时间步长可足够小（例如0.01）。给出t=20,50两个时刻的数值解，与精确解比较（画图），并给出数值解的L2模误差。,34,Copyright by Li Xinliang,提示：1.如不使用优化技术，则格式构造方法简单，Taylor展开后解代数方程组即可。2.建议尝试使用优化技术,例：假设格式形式如下,如果要求其有5阶精度，则通过Taylor展开可得到6个方程，6个系数可直接解出。我们要求其有4阶精度（当然3阶，2阶也可），于是Taylor展开只能提供5个方程。6个未知数（a1-a6），5个方程；有1个自由参数。调整这个自由参数，使得kr,ki曲线最为理想。如何调整？1)可以人工调整，观

22、察kr,ki曲线，选取满意的。2）可自动调整，设立一个优化目标函数。例如 调整自由参数，使得该目标函数取最大值。思路：牺牲精度，提高分辨率,35,Copyright by Li Xinliang,作业题2（选作题）：构造更高分辨率的GVC格式,对于空间导数，构造出一种不超过6点的GVC格式。要求：a.精度不限；b.网格基架点数不超过6个；c.求解模型方程,计算结果间断尽量保持“锐利”;计算结果振荡尽量小。,振荡的定量判据：总变差（Total Variation）：间断“锐利”的定量判据：间断区内的点数？（自行设计）,给出差分格式的表达式、色散/耗散分析（ki,kr曲线）；给出模型方程t=0.2

23、的结果（空间100个网格点，计算域0,1，时间推进可采用3阶Runge-Kutta方法）；与精确解及NND2a进行比较（画在同一张图上）,建议：利用优化方法,36,Copyright by Li Xinliang,作业题3：数值求解Sod 激波管问题,计算其数值解，画出t=0.14时刻密度、速度及压力的分布；并与精确解进行比较（要求画在一张图上）。要求：1）空间网格数100，时间推进格式选用3阶Runge-Kutta,时间步长自选。2）可选用逐点分裂，也可选用特征分裂。3）建议采用本讲作业题2（或作业题1）自行构造的差分格式计算。（作业题2是激波捕捉格式，效果应当会好些）。如果作业题1和作业题2遇到困难，也可采用现有的差分格式。,37,Copyright by Li Xinliang,