计算流体力学-不可压缩N-S方程的求解.课件.ppt

1、计算流体力学讲义计算流体力学讲义 第十二讲第十二讲 不可压缩不可压缩Navier-Stokes方程的求解方程的求解李新亮李新亮 ；力学所主楼；力学所主楼219； 82543801 知识点：知识点： 1讲义、课件上传至讲义、课件上传至 (流体中文网）流体中文网） - “流体论坛流体论坛” -“ CFD基础理论基础理论 ”讲课录像及讲义上传至网盘讲课录像及讲义上传至网盘 http:/cid- by Li Xinliang拟压缩性方法拟压缩性方法求解压力求解压力Poisson方程法方程法涡流函数法涡流函数法Simple方法方法Copyright by Li Xinliang2知识回顾知识回顾一、一、

2、 代数方程组的求解代数方程组的求解nnnnnnnnaaaaaaaaaaaaaaaaA321333323122322211131211bAx 直直接接法法nnnnnaaaaaaaaaaA0000000333223221131211Gauss消元法nnnnnnnnnnnnnnnnuuuuuuuuuullllllaaaaaaaaaaaaaaaa3332232211312113213231213213333231223222111312111111LU分分解法解法追赶法：追赶法：jjjjjjjdxcxbxa11jjjjBxAx1Copyright by Li Xinliang3迭代法迭代法),(),(

3、2222yxguyxfyuxu2,1,1, 1, 14jijijijijijifuuuuuJacobi迭代迭代n+1nnnnGauss-Seidel迭代迭代n+1nn+1n+13n+1n+1nn+1n+1n+1nnnLU-ADILU-SGS知识回顾知识回顾Copyright by Li Xinliang4知识回顾知识回顾二、二、 网格生成网格生成)(22xyy )(11xyy ABCDEFABCDEF物理空间计算空间1. 代数网格生成法代数网格生成法2. 解椭圆型方程网格生成法解椭圆型方程网格生成法Copyright by Li Xinliang512. 1 不可压缩不可压缩Navier-St

4、okes方程的特点方程的特点VVVVV2Re10pt密度为常数的不可压缩密度为常数的不可压缩Navier-Stokes方程组：方程组：特点：特点： 动量方程与能量方程解耦动量方程与能量方程解耦 压力属于压力属于约束变量约束变量而不是而不是发展变量发展变量温度对密度的影温度对密度的影响可忽略不计响可忽略不计压力不能时间推进压力不能时间推进求解求解可压缩可压缩不可压缩不可压缩易处易处压力可推进求解，易于使用压力可推进求解，易于使用显格式显格式压力方程具有椭圆性，无法推进压力方程具有椭圆性，无法推进求解。压力方程收敛性差求解。压力方程收敛性差难处难处可能出现间断可能出现间断不会出现间断不会出现间断研

5、究重点研究重点激波捕捉激波捕捉压力处理压力处理)(VVVVCopyright by Li Xinliang6概念澄清：概念澄清： 压力压力 动力学压力动力学压力及及热力学压力热力学压力动力学压力动力学压力 应力的中各向同性部分应力的中各向同性部分ijijijpP211p22p12p22pnPpn连续介质微元体的受力平衡：连续介质微元体的受力平衡： 应力的概念应力的概念nppnpp静止流体或无粘流体中力的静止流体或无粘流体中力的平衡平衡 动力学压力的概念动力学压力的概念热力学压力热力学压力 由分子动力学性质决定由分子动力学性质决定 状态方程状态方程热力学压力：热力学压力：分子对固壁的碰撞，分子对

6、固壁的碰撞，产生压力产生压力RTp完全气体：动力学压力热力学压力可压缩可压缩N-S方程：方程： 动力学与热力学耦合；动力学压力动力学与热力学耦合；动力学压力= 热力学压力热力学压力不可压缩不可压缩N-S方程：方程： 动力学与热力学解耦动力学与热力学解耦 由不可压缩条件确定压力由不可压缩条件确定压力 （纯动力学概念）（纯动力学概念） 1） 压力的处理原则压力的处理原则Copyright by Li Xinliang7 奇偶失联与交错网格奇偶失联与交错网格uypyvvxvutvuxpyuvxuutuyvxu22Re1Re10压力项，通压力项，通常采用中心常采用中心差分离散差分离散yppypxppx

7、pjijijijijiji2,21,1, 1, 1,极端情况：极端情况： 棋盘式压力场棋盘式压力场 高压高压 低压低压02,021,1, 1, 1yppxppjijijiji特点：特点： 高压高压-低压点间隔分布低压点间隔分布采用中心差分格式计算出：采用中心差分格式计算出：0ypxp流场竟然流场竟然“保持稳定保持稳定”“奇偶失联奇偶失联”Copyright by Li Xinliang8常用措施：常用措施： 交错网格交错网格 压力压力p 速度速度 uuypyvvxvutvuxpyuvxuutuyvxu22Re1Re10 速度速度v交错网格示意图交错网格示意图jip,jiu,2/12/1, ji

8、vuxpyuvxuutu2Re1在u的网格点上离散jip,jiu,2/12/1, jivxppxpjijiji, 1,2/1uypyvvxvutv2Re1在v的网格点上离散xppypjijiji,1,2/1,jip, 1注：注： 对流项通常采用迎风格式离散对流项通常采用迎风格式离散xuuxuuxuu 后差 前差2uuuCopyright by Li Xinliang92) 对流项的处理原则对流项的处理原则 VVVVV2Re10ptVVVVVVVV)()()()()(yvxuuyyvxuuyuvxuu关系式1：)()()(yvxuvyvvxvuyvvxvu关系式2：VVV221VVVV2Re10

9、pVt兰姆兰姆-葛罗米柯等式葛罗米柯等式2/2Vpp总压表达式表达式优点优点不足不足普通型普通型简单，易于迎风简单，易于迎风有混淆误差有混淆误差守恒型守恒型简单，守恒简单，守恒有混淆误差有混淆误差旋度型旋度型混淆误差小混淆误差小计算量大计算量大螺旋螺旋-对称型对称型混淆误差小混淆误差小计算量略大计算量略大普通型守恒型VV )(VVV2/)(VVVVVCopyright by Li Xinliang1012. 1 人工压缩性方法（求解定常方程）人工压缩性方法（求解定常方程）VVVVV2Re10ptVVVVV2Re10pttp人工压缩性因子人工压缩性因子达到定常态达到定常态0 V0流动压缩时流动压

10、缩时 （ ），压力升高），压力升高流动膨胀时流动膨胀时 （ ），压力降低），压力降低0 V0 V0Vtp0 V0 V 增大增大 可令压力收敛加快，可令压力收敛加快，但会增加方程的刚性（降低时但会增加方程的刚性（降低时间步长）间步长）。人工压缩性因子 相当于2ctpctppts202VctpCopyright by Li Xinliang11对于定常问题，需要迭代到收敛对于定常问题，需要迭代到收敛VVVVV2Re10pttp)1,max(111nnnnnnppvvuu对于非定常问题，需要内迭代对于非定常问题，需要内迭代 （效率较低）（效率较低）VVVVV2Re10ptVVVVV2Re10pttp

11、nknnkkkkkpttppVVVVVV211Re10Step 1 : 得到得到n 时间步的值时间步的值nVStep 2: 进行如下内迭代直至收敛进行如下内迭代直至收敛Step 3: 收敛后的收敛后的V即为即为1nV内迭代收敛慢，效率较低；内迭代收敛慢，效率较低；通常不使用人工压缩方法解非通常不使用人工压缩方法解非定常问题。定常问题。Copyright by Li Xinliang1212.2 求解压力求解压力Poisson方法方法 （投影法）（投影法）1） 压力的控制方程压力的控制方程VVVVV2Re10pt对动量方程求散度对动量方程求散度)(2VV pVVVVVV22Re1)(ptpPoi

12、sson方程方程压力的控制方程压力的控制方程无法时间推进无法时间推进需联立求解，通常采用需联立求解，通常采用时间分裂法时间分裂法Copyright by Li Xinliang132） 投影法投影法求解微分型压力求解微分型压力Poisson方程方程0)Re1(02ptVVVVV原理：原理： 将时间推进分成三个子步，将时间推进分成三个子步， 中间步解出压力中间步解出压力可时间推进可时间推进不能时间不能时间推进推进Step 1: 预算步预算步0)Re1(2*nntVVVVVStep 2: 压力修正步压力修正步001*1nnptVVV*21Vtp求解，得到压力求解，得到压力pStep 3: 最终步最

13、终步0*1ptnVV得到得到n+1时刻的时刻的V以以 1阶精度时间推进方法为例，实际上可采用更高阶精度时间推进方法：阶精度时间推进方法为例，实际上可采用更高阶精度时间推进方法： Karniadakis GE, Israeli M, Orszag SA. 1991 High-order splitting methods for the incompressible Navier-stokes equations. J. Comp. Phys. 97:414-443. Copyright by Li Xinliang143） 投影法投影法求解离散型压力求解离散型压力Poisson方程方程压力修正

14、步：压力修正步： 将将离散离散的动量方程带入的动量方程带入离散离散的连续性方程，得到的连续性方程，得到离散离散的压力方程的压力方程001*1nnptVVV)3(0/ )(/ )()2(0)(/) 1 (0)(/12/1,12/1,1,2/11,2/11,1,*2/1,12/1, 1*,2/11,2/1yvvxuuppytvvppxtuunjinjinjinjijijijinjijijijinji交错网格上离散将将 （1） （2） 两式两式 （离散的动量方程）（离散的动量方程）带入带入 （3） 式式 （离散的连续性方程）可得到关于压力（离散的连续性方程）可得到关于压力p的方程的方程 （离散的压力

15、（离散的压力Poisson方程）；方程）；该方法可保证该方法可保证 （3） 式严格满足，因而相容性比方法式严格满足，因而相容性比方法2）更好）更好Copyright by Li Xinliang1512.3 涡量涡量-流函数方法流函数方法 （二维问题）（二维问题）)3(Re1)2(Re1) 1 (022uypyvvxvutvuxpyuvxuutuyvxu)3()2(xyxvyu)4(Re12yvxut引入流函数引入流函数xvyu,)5(2（4） （5） 式即涡量式即涡量-流函数的控制方程流函数的控制方程计算结束后，如果需要计算压力，则求解如下方程计算结束后，如果需要计算压力，则求解如下方程22

16、222yvxvyuxup)3()2(yxCopyright by Li Xinliang16驱动方腔流动驱动方腔流动例：例： 求解驱动方腔流动求解驱动方腔流动问题描述：问题描述： 如图示边长为如图示边长为L的方腔，上表面流体以常速度的方腔，上表面流体以常速度U运动，求解里面的流场（假设流动定常）。运动，求解里面的流场（假设流动定常）。以涡量以涡量-流函数法为例：流函数法为例：2Re1yvxut400ReULUL21） 离散化离散化jixxjiuuxu,)(对流项：对流项：迎风差分迎风差分2uuu迎风差分，建议采用高阶的迎风差分，建议采用高阶的粘性项：采粘性项：采用中心差分用中心差分2, 1,

17、1222xxjijiji也可采用更高阶的，可借也可采用更高阶的，可借助求差分系数的小程序助求差分系数的小程序时间推进：时间推进： 可采用显格式可采用显格式ttnn1Copyright by Li Xinliang172采用中心差分离散：采用中心差分离散：jijijijijijijiyx,2,1,1,2, 1, 122可采用可采用Jocabi，Gauss-Seidel等方法迭代等方法迭代提示：提示： 时间推进过程中的中间步无需迭代至收敛，时间推进过程中的中间步无需迭代至收敛， 最终（最后一个时最终（最后一个时间步）收敛即可。间步）收敛即可。2） 边界条件边界条件UL速度边界条件：速度边界条件：

18、上壁面上壁面u=1,v=0; 其他壁面其他壁面u=v=0;流函数的边界条件：流函数的边界条件： 0边界是一条流线，边界是一条流线， 流线是流函数的等值线流线是流函数的等值线涡量的边界条件：涡量的边界条件： 由速度给出由速度给出niniiyyuxuuxyu)(1)(12,1 ,2,00可用更高阶的格式可用更高阶的格式Copyright by Li Xinliang1812.4 SIMPLE方法方法基本思想：基本思想： 与（离散型）投影法类似，与（离散型）投影法类似， 但速度推进是隐式的；但速度推进是隐式的；0)Re1(02ptVVVVV0)Re1(*2*ptnVVVVV1） 已知预估压力已知预估

19、压力 计算速度计算速度*p采用隐式采用隐式离散离散) 1 ()(.)(*,*, 1*,2/32*,2/11*,2/1jijijijijippuauau)2()(.)(*,*1,*2/3,2*2/1,1*2/1,jijijijijippvbvbv2） 压力及速度修正压力及速度修正*vvvuuuppp)(.)(, 1,2/32,2/11,2/1jijijijijippuauau)(.)(,1,2/3,22/1,12/1,jijijijijippvbvbv重要简化重要简化)(, 1,2/1jijijippu)(,1,2/1,jijijippv类似类似“人工压缩性方法人工压缩性方法”修正方程对角化修正

20、方程对角化 （显式化）（显式化）已知已知联立求解Copyright by Li Xinliang19带入离散的连续性方程带入离散的连续性方程：0/ )(/ )(12/1,12/1,1,2/11,2/1yvvxuunjinjinjinji)(, 1*,2/11,2/1jijijinjippuu)(,1,*2/1,12/1,jijijinjippvv得到离散的压力得到离散的压力Poisson方程：方程：)5(1,41,3, 12, 11,cpcpcpcpcpjijijijiji求解后，得到压力修正值：求解后，得到压力修正值：jip,（4）带入（带入（4）时得到）时得到n+1时刻的速度时刻的速度具体

21、步骤：具体步骤： 1） 已知已知n时刻的速度压力时刻的速度压力 2） 预估压力预估压力 （可取为（可取为n时刻的压力）时刻的压力） 3） 带入（带入（1）（）（2）式，解出）式，解出 （隐格式隐格式，需迭代求解），需迭代求解） 4） 求解压力的修正方程求解压力的修正方程 （5）得到修正压力）得到修正压力 5） 带入（带入（4）式，得到）式，得到n+1时刻的速度及压力时刻的速度及压力 6） 推进求解直到给定时刻（或收敛）推进求解直到给定时刻（或收敛）*p*,vu如该步改用如该步改用显格式，则显格式，则为（离散型）为（离散型）投影法投影法*1pppn提示：提示： 对于定常问题，内迭代无需收敛，最终

22、时刻收敛即可对于定常问题，内迭代无需收敛，最终时刻收敛即可习题习题 12.1 求解方腔问题求解方腔问题问题描述：问题描述： 如图示边长为如图示边长为L的方腔，上表面流体以常速度的方腔，上表面流体以常速度U运运动，求解里面的流场（假设流动定常）。动，求解里面的流场（假设流动定常）。 考虑考虑 三种情况三种情况1000,400,100ReULUL要求：要求： 数值方法不限数值方法不限 （人工压缩性方法、投影法、涡量（人工压缩性方法、投影法、涡量-流函数方法及流函数方法及SIMPLE方方法均可）；法均可）； 空间离散采用差分法，建议采用较高阶精度的方法。空间离散采用差分法，建议采用较高阶精度的方法。 绘制出定常解的流线图。绘制出定常解的流线图。 请详细写明方程及公式的推导过程及计算流程，切勿只上交计算结果。请详细写明方程及公式的推导过程及计算流程，切勿只上交计算结果。