第13章_计算流体力学CFD(5)

1、Lax-Wendroff方法方法是一种显式有限差是一种显式有限差分方法，适合于推分方法，适合于推进求解。进求解。二维时间推进网格二维时间推进网格Lax-Wendroff方法方法在时间和空间上都在时间和空间上都具有二阶精度。具有二阶精度。二维时间推进网格二维时间推进网格非定常二维无粘流（欧拉方程非守恒形式）：非定常二维无粘流（欧拉方程非守恒形式）：Lax-Wendroff显式推进求解显式推进求解 (沿时间方向进行泰勒级沿时间方向进行泰勒级数展开数展开)：空间导数采用中心差分：空间导数采用中心差分：()求对时间求对时间t的二阶导数：的二阶导数：Lax-Wendroff显式推进求解显式推进求解 ：M

2、acCormack方法在时间和空间上都方法在时间和空间上都具有二阶精度。具有二阶精度。MacCormack方法是一种显式有限差分方法是一种显式有限差分方法，适合于推进求解。方法，适合于推进求解。MacCormack方法比方法比Lax-Wendroff方方法应用起来更简单。法应用起来更简单。校正步校正步预估步预估步预估步：空间导数用向前差分计算。预估步：空间导数用向前差分计算。预估步：空间导数用向前差分计算。预估步：空间导数用向前差分计算。预估值：预估值：校正步：空间导数用向后差分计算。校正步：空间导数用向后差分计算。在在MacCormack方法中，预估步用向前差分，方法中，预估步用向前差分，校

3、正步用向后差分；也可以预估步用向后差分，校正步用向后差分；也可以预估步用向后差分，校正步用向前差分。或者在时间推进解法的相校正步用向前差分。或者在时间推进解法的相继两个时间步中轮流使用这两种办法。继两个时间步中轮流使用这两种办法。粘性流动的控制方程是粘性流动的控制方程是N-S方程。方程。对定常流动，对定常流动，N-S方程的数学性质更多地表现为方程的数学性质更多地表现为椭圆型的，不能采用椭圆型的，不能采用Lax-Wendroff方法和方法和MacCormack方法求解。方法求解。对非定常流动，可以采用对非定常流动，可以采用Lax-Wendroff方法或方法或MacCormack方法求解方法求解N

4、-S方程。方程。22vevuU可以采用可以采用Lax-Wendroff方法或方法或MacCormack方法求解方法求解U的分的分量在各时间步的值。量在各时间步的值。非定常守恒形式欧拉方程（二维）：非定常守恒形式欧拉方程（二维）：定常守恒型二维欧拉方程：定常守恒型二维欧拉方程：对于亚声速流动，上述对于亚声速流动，上述方程是椭圆型的，所有方程是椭圆型的，所有空间推进方法都不适用，空间推进方法都不适用，MacCormack方法也不方法也不适用。适用。对于超声速流动，上述方对于超声速流动，上述方程是双曲型的，空间推进程是双曲型的，空间推进方法适用，方法适用，MacCormack方法也适用。方法也适用。

5、定常守恒型二维欧拉方程：定常守恒型二维欧拉方程：MacCormack方法：方法：定常守恒型二维欧拉方程：定常守恒型二维欧拉方程：预测步预测步：（向前差分）：（向前差分）预估值：预估值：预估值：预估值：校正步校正步：（向后差分）：（向后差分）松弛法特别适合于求解椭圆型偏微分方程，松弛法特别适合于求解椭圆型偏微分方程，常被用来求解无粘亚声速的低速流动。常被用来求解无粘亚声速的低速流动。考虑无粘不可压流体的二维无旋流动，控考虑无粘不可压流体的二维无旋流动，控制方程为制方程为Laplace方程：方程：松弛法是一种迭代法松弛法是一种迭代法上标上标n和和n+1表示迭代次数表示迭代次数松弛法是一种迭代法松弛

6、法是一种迭代法松弛法是一种迭代法松弛法是一种迭代法松弛法是一种迭代法松弛法是一种迭代法从左至右扫描从左至右扫描松弛法是一种迭代法松弛法是一种迭代法当所有网格点处的当所有网格点处的 都小于一个预定的值时，迭代都小于一个预定的值时，迭代收敛。收敛。njinji,1,运用逐次松弛法可加运用逐次松弛法可加快收敛的过程。快收敛的过程。从左至右扫描从左至右扫描从下至上扫描从下至上扫描运用逐次松弛法可加运用逐次松弛法可加快收敛的过程。快收敛的过程。 是松弛因子，如果是松弛因子，如果 1，叫做逐次超松弛法；，叫做逐次超松弛法；如果如果 1，叫做逐次低松弛法。，叫做逐次低松弛法。运用逐次松弛法可加运用逐次松弛法

7、可加快收敛的过程。快收敛的过程。选取合适的选取合适的 值，可以减少迭代次数，从而减少计算值，可以减少迭代次数，从而减少计算时间。在某些问题中，迭代次数可减少到原来的时间。在某些问题中，迭代次数可减少到原来的1/30一维波动方程：一维波动方程：差分方程：差分方程：截断误差：截断误差：差分方程：差分方程：泰勒级数展开：泰勒级数展开：差分方程：差分方程：将泰勒级数展开代入差分方程得：将泰勒级数展开代入差分方程得：差分方程：差分方程：将泰勒级数展开代入差分方程得：将泰勒级数展开代入差分方程得：差分方程：差分方程：等号右边将对等号右边将对t的偏导数转化为对的偏导数转化为对x的偏导数得：的偏导数得：差分方

8、程：差分方程：偏微分方程（修正方程）：偏微分方程（修正方程）：一维波动方程（偏微分方程）：一维波动方程（偏微分方程）：差分方程：差分方程：一维波动方程（偏微分方程）：一维波动方程（偏微分方程）：差分方程的精确解是上述一维波动方程的数值解（含误差）差分方程的精确解是上述一维波动方程的数值解（含误差）差分方程：差分方程：差分方程的精确解是上述修正方程的精确解（不含误差）差分方程的精确解是上述修正方程的精确解（不含误差）偏微分方程（修正方程）：偏微分方程（修正方程）：差分方程：差分方程：偏微分方程（修正方程）：偏微分方程（修正方程）：一维波动方程（偏微分方程）：一维波动方程（偏微分方程）：修正方程等

9、号右端的项是截断误差，如果截断误差的主项修正方程等号右端的项是截断误差，如果截断误差的主项是偶数阶导数，数值解将主要表现出耗散行为；如果主项是偶数阶导数，数值解将主要表现出耗散行为；如果主项是奇数阶导数，数值解将主要表现出色散行为。是奇数阶导数，数值解将主要表现出色散行为。偏微分方程（修正方程）：偏微分方程（修正方程）：等号右端的偶数阶导数项起数值耗散的作用，奇数阶导数等号右端的偶数阶导数项起数值耗散的作用，奇数阶导数项起数值色散的作用。项起数值色散的作用。偏微分方程（修正方程）：偏微分方程（修正方程）：数值耗散的作用很象物理粘性，二阶导数项前的系数被称数值耗散的作用很象物理粘性，二阶导数项前

10、的系数被称为人工粘性。为人工粘性。偏微分方程（修正方程）：偏微分方程（修正方程）：数值耗散的影响会将波抹平数值耗散的影响会将波抹平色散导致波的不同相位在传播中产生畸变，色散导致波的不同相位在传播中产生畸变，表现为波前和波后出现振荡。表现为波前和波后出现振荡。尽管人工粘性降低了解的精度，但通常有助于提高解的稳尽管人工粘性降低了解的精度，但通常有助于提高解的稳定性。定性。偏微分方程（修正方程）：偏微分方程（修正方程）：考虑二维热传导方程：考虑二维热传导方程：等号右端有五个未知量，不能得到三对角方程组，不能采等号右端有五个未知量，不能得到三对角方程组，不能采用托马斯算法（追赶法）求解。用托马斯算法（

11、追赶法）求解。采用采用Crank-Nicolson方法（隐式）：方法（隐式）：考虑二维热传导方程：考虑二维热传导方程：第一步：第一步：时间步长为时间步长为 ，空间导数采用中心差分，只对，空间导数采用中心差分，只对x的的导数采用隐式处理。导数采用隐式处理。第一步：第一步：简化为三对角形式简化为三对角形式第一步：第一步：对每一个固定的对每一个固定的j，对所有，对所有的的i联立形成方程组。联立形成方程组。对不同的对不同的j，重复上述过程。，重复上述过程。考虑二维热传导方程：考虑二维热传导方程：第二步：第二步：时间步长为时间步长为 ，空间导数采用中心差分，只对，空间导数采用中心差分，只对y的的导数采用

12、隐式处理。导数采用隐式处理。第二步：第二步：简化为三对角形式简化为三对角形式第二步：第二步：对每一个固定的对每一个固定的i，对所有，对所有的的j联立形成方程组。联立形成方程组。对不同的对不同的i，重复上述过程。，重复上述过程。两步结束之后，两步结束之后，T在时间方向上推进了一个时间步长在时间方向上推进了一个时间步长 t.考虑二维热传导方程：考虑二维热传导方程：推进过程只涉及三对角方程组。推进过程只涉及三对角方程组。第一步，差分方程的第一步，差分方程的x方向是隐式的。方向是隐式的。考虑二维热传导方程：考虑二维热传导方程：所以这种方法叫交替方向隐式方法所以这种方法叫交替方向隐式方法(Alterna

13、ting Direction Implicit, ADI)第二步，差分方程的第二步，差分方程的y方向是隐式的。方向是隐式的。考虑二维热传导方程：考虑二维热传导方程：ADI格式对格式对t,x,y都是二阶精度的都是二阶精度的截断误差为：截断误差为：不可压无粘流动受椭圆型偏微分方程控制（不可压欧不可压无粘流动受椭圆型偏微分方程控制（不可压欧拉方程），松弛法是求解椭圆型问题经典的数值方法，拉方程），松弛法是求解椭圆型问题经典的数值方法，本质上是一个迭代过程。本质上是一个迭代过程。不可压粘性流动的控制方程是不可压的不可压粘性流动的控制方程是不可压的NS方程，这方程，这个方程具有椭圆型和抛物型的混合特性，

14、松弛法不是个方程具有椭圆型和抛物型的混合特性，松弛法不是特别适用。特别适用。压力修正法也是一种迭代过程，在不可压压力修正法也是一种迭代过程，在不可压NS方程的方程的数值求解中得到了广泛的应用。数值求解中得到了广泛的应用。假设假设 =常数，常数， =常数，可压缩常数，可压缩NS方程转化为不可方程转化为不可压压NS方程：方程：上述四个方程封闭，含上述四个方程封闭，含 四个未知数。四个未知数。二维不可压流体的连续二维不可压流体的连续性方程为：性方程为：中心差分格式为：中心差分格式为：右上角是右上角是u的值，的值，左下角是左下角是v的值的值速度会出现右图的速度会出现右图的棋盘式分布棋盘式分布右上角是右

15、上角是u的值，的值，左下角是左下角是v的值的值可压流动中不会发生右可压流动中不会发生右图的问题，因为连续性图的问题，因为连续性方程中包含了密度对时方程中包含了密度对时间和空间的变化。间和空间的变化。0Vt在可压缩流动中，右图在可压缩流动中，右图速度的棋盘分布经过一速度的棋盘分布经过一个时间步就会被抹平。个时间步就会被抹平。二维不可压流体压力梯二维不可压流体压力梯度采用中心差分：度采用中心差分：棋盘式的离散压力分棋盘式的离散压力分布布压力会出现右图的压力会出现右图的棋盘式分布棋盘式分布在交错网格上使用中心在交错网格上使用中心差分就不会出现速度和差分就不会出现速度和压力的棋盘式分布问题。压力的棋盘

16、式分布问题。交错网格交错网格在在(i-1,j), (i,j), (i+1,j), (i,j+1),(i,j-1)等图中的实等图中的实心原点上计算压力心原点上计算压力交错网格交错网格在在(i-1/2,j), (i+1/2,j)等图等图中的空心原点上计算中的空心原点上计算u交错网格交错网格在在(i,j-1/2), (i,j+1/2)等图等图中的空心原点上计算中的空心原点上计算v连续性方程在网格点连续性方程在网格点(i,j)的中心差分表达式为：的中心差分表达式为：交错网格交错网格压力修正法本质上是一种迭代法，思路如下：压力修正法本质上是一种迭代法，思路如下：1）迭代开始时，先给定压力的初始近似）迭代

17、开始时，先给定压力的初始近似p*2）用）用p*的值从动量方程中求解的值从动量方程中求解u,v,w,得到与得到与p*有关有关的的u*,v*,w*压力修正法本质上是一种迭代法，思路如下：压力修正法本质上是一种迭代法，思路如下：修正后的压力为修正后的压力为3）将）将u*,v*,w*代入连续性方程，它们不一定满足连代入连续性方程，它们不一定满足连续性方程。用连续性方程构造压力的修正量续性方程。用连续性方程构造压力的修正量 ,加加到到p*上，使速度场满足连续性方程。上，使速度场满足连续性方程。修正后的速度为修正后的速度为速度修正量速度修正量 可以从可以从 得到。得到。压力修正法本质上是一种迭代法，思路如

18、下：压力修正法本质上是一种迭代法，思路如下：4） 用步骤用步骤3)中修正后的压力做为新的中修正后的压力做为新的p*，回到步，回到步骤骤2)。重复这个过程，直到速度场满足连续性方程。重复这个过程，直到速度场满足连续性方程为止。为止。这样就得到修正好了的流场。这样就得到修正好了的流场。压力修正公式为：压力修正公式为：压力修正公式为：压力修正公式为：上述压力修正公式具有椭圆型的性质，可以用松弛上述压力修正公式具有椭圆型的性质，可以用松弛法数值求解。法数值求解。在不可压流场中，压力的扰动将会传遍整个流场，在不可压流场中，压力的扰动将会传遍整个流场，这与上述方程的椭圆型性质相吻合。这与上述方程的椭圆型性

19、质相吻合。压力修正公式为：压力修正公式为：压力修正公式是压力修正压力修正公式是压力修正 的泊松方程的中心差的泊松方程的中心差分表达式。分表达式。上述泊松方程中的二阶偏导数用中心差分替代。上述泊松方程中的二阶偏导数用中心差分替代。tdQ/式中：式中：压力修正压力修正 的泊松方程的泊松方程(为椭圆型为椭圆型)：tdQ/d相当于一个质量源项。相当于一个质量源项。SIMPLE是是Semi-implicit method for pressure-linked equation (压力耦合方程的半隐式算法压力耦合方程的半隐式算法)的缩写。的缩写。SIMPLE算法的步算法的步骤如下：骤如下：1）在右图所示

20、的交）在右图所示的交错网格上分别给出错网格上分别给出*,np*,nu*nvSIMPLE算法的步算法的步骤如下：骤如下：2）求出）求出1*,nu1*nv采用动量方程求解。采用动量方程求解。2） 的求法：的求法：1*nuX方向的动量方程：方向的动量方程：2） 的求法：的求法：1*nu在在a点：点：在在b点：点：X方向的动量方程：方向的动量方程：差分方程：差分方程：差分方程：差分方程：X方向的动量方程：方向的动量方程：2） 的求法：的求法：1*nu2） 的求法：的求法：1*nu2） 的求法：的求法：1*nvY方向的动量方程：方向的动量方程：2） 的求法：的求法：1*nv在在c点：点：在在d点：点：Y方向的动量方程：方向的动量方程：差分方程：差分方程：2） 的求法：的求法：1*nv3）将）将 和和 代入压力修正公式，在所有内代入压力修正公式，在所有内部网格点上求解部网格点上求解 1*nu1