《湍流及转捩》PPT课件.pptVIP

计算流体力学讲义 第十三讲 湍流与转捩 （1） 李新亮 ；力学所主楼219； 82543801,知识点：,1,讲义、课件上传至 (流体中文网） - “流体论坛” -“ CFD基础理论 ” 讲课录像及讲义上传至网盘 /browse.aspx/.Public,Copyright by Li Xinliang,线性稳定性理论 转捩的预测方法 壁湍流转捩的涡动力学机制,Copyright by Li Xinliang,2, 13.1 线性稳定性理论,一、 稳定性基本概念,常识：流体中的不稳定性,K-H不稳定性,A. K-H (Kelvin-Helmholtz）不稳定性 自由剪切流的无粘不稳定性,混合层 K-H不稳定性,K-H不稳定性的关键： 速度剖面有拐点,Lee-Lin: 速度剖面的拐点是无粘不稳定性的必要条件,流体不禁搓，一搓搓出涡,已知某运动状态； 在此基础上施加微小扰动； 如扰动随时间（或空间）衰减，则称系统稳定，否则为不稳定,Copyright by Li Xinliang,3,自然界中 K-H不稳定性图片,智利塞尔扣克岛的卡门涡街,澳大利亚Duval山上空的云,KelvinHelmholtz instability clouds in San Francisco,佛兰格尔岛周围的卡门涡街,高速流,低速流,自由剪切层受到扰动界面变形后的情况 K-H不稳定性的产生机理,受阻减速，压力升高，产生高压区,高压导致变形加剧,Copyright by Li Xinliang,4,B. T-S (Tollmien-Schlichting) 不稳定性不可压 壁面剪切流的粘性不稳定性,Mack 不稳定性 超声速壁面剪切流的不稳定性,不可压边界层速度剖面 （Blasius解） 无拐点,可压缩情况 Mach数足够高时会出现广义拐点 出现无粘不稳定性,不可压缩无粘不稳定性 需存在拐点 可压缩无粘不稳定性需存在广义拐点,Mach 6 钝锥（1攻角） 不同子午面 的分布,超音速平板边界层的不稳定波,第1模态（T-S波）,第2模态 （Mack模态）,Copyright by Li Xinliang,5,激波,密度界面,R-M (Richtmyer-Meshkov)不稳定性 激波与密度界面作用的斜压效应,惯性约束聚变（ICF）示意图,小知识 涡的产生机制： 粘性、 斜压、有旋的外力,激波,密度界面,斜压项,Copyright by Li Xinliang,6,D. R-T (Reyleigh-Taylor)不稳定性 重力带来的不稳定性,R-T (Reyleigh-Taylor)不稳定性,重介质,轻介质,Copyright by Li Xinliang,7,E Barnard热对流不稳定性,其他学科的不稳定性：,Euler压杆的不稳定性,Barnard 热对流的胞格结构,板壳的不稳定性,Copyright by Li Xinliang,8,二、 稳定性问题的常用数学方法 线性稳定性分析,Step 1: 得到线性化的扰动方程,控制方程为：,已知其具有解,最好是精确解，也可用高精度的数值解,令：,舍弃高阶小量，得到线性化的扰动方程,（1）,例如： 平板的Blasius解，槽道的Poiseuille 解,线性方程,Copyright by Li Xinliang,9,Step 2: 求解 的特征值问题,什么条件下具有非零解，非零解如何？,通常假设在某些方向具有周期性，转化为一维问题,数值方法： 将 （1） 离散代数方程 何时有非零解， 非零解如何？ 特征值问题,什么条件下有非零解？,特征值问题计算量巨大，目前通常只能处理一维问题,Copyright by Li Xinliang,10,三、 稳定性问题示例 不可压缩槽道流动的线性稳定性（LST）理论 (以二维为例),Step 1: 获得线性化扰动方程,令：,Poiseuille解：,(2),代入方程（2），并舍去高阶小量得到线性化的扰动方程,（3）,1） 控制方程及边界条件,Copyright by Li Xinliang,11,研究扰动发展的空间模式和时间模式,扰动源,空间模式： 任一点的扰动具有时间周期性 符合物理条件,假设扰动具有如下形式：,沿流向及时间方向具有波动特性 称为Tolmien-Schlichting（T-S）波,任意扰动可分解为正弦波的叠加 线性系统各成分无相互作用 可独立研究,为实数,为复数,扰动波的振幅沿流向指数变化,空间增长率,时间模式： 扰动具有流向的周期性 假设一窗口沿流向运动，研究窗口内扰动的演化,为实数,为复数,扰动波的振幅虽时间变化,时间增长率,Copyright by Li Xinliang,12,以时间模式为例：,（4）,（5）,（6）,线性偏微方程（3）转化成为含参数的线性常微方程组（4）-（6）,谱方法的常规做法,通过消元法，转化为更高阶的常微方程 （不是必须的）,常用做法，通常还可以反向为之： 高阶方程转化为低阶方程组,消去,Orr-Sommerfeld（O-S) 方程,其中：,最终，控制方程为O-S方程：,Copyright by Li Xinliang,13,边界条件：,y=1 （固壁）:,y=0 （中心线，对称）：,可以取计算域-1,1,使用固壁边界条件； 也可以取计算域-1,0,使用固壁及对称边界条件,流函数形式的O-S方程,引入流函数，使得：,计算出 后，利用公式,计算其他两个量,则：,令：,常数倍,满足的方程及边界条件与 完全相同。,如果 恒大于（或恒小于0），则必有,Copyright by Li Xinliang,14,小知识： 关于O-S方程,1） O-S方程适用于不可压平行流的稳定性问题 （不仅槽道流） 2） 准平行流 （流线沿x方向接近平行）也可使用（例如边界层流动） 3）如果舍去粘性（左端）项，则方程称为Rayleigh方程,Rayleigh拐点定理： Rayleigh方程存在不稳定解的必要条件是速度型存在拐点。即存在某点 使得,若存在无粘不稳定性，该项必有0点。,分部积分，并取虚部，得：,不存在非稳定解,Copyright by Li Xinliang,15,2） O-S方程的解法,数学表述 奇性（特征值）问题： 参数 为何值时，方程有非零解？ 非零解如何？,时间发展槽道湍流： (通常) 给定Re及a，问 w取何值时，O-S方程有非零解？,增长率,求解步骤： 1） 将O-S方程离散，得到线性代数方程组 离散方法： 差分法、有限元法、谱方法、打靶法 2） 求w，使得该方程有非零解（奇性或特征值问题）.,求出w,局部法：只求出一个w 全局法：计算出全部的w,试计算 的时间发展槽流中 （即波长2p）T-S波的频率及增长率,Copyright by Li Xinliang,16,四: 例题,Step 1: 离散 （差分法）,一维问题,网格： 均匀网格 （简单，但需要较多网格点 N 300） 非均匀网格,差分离散：,2阶格式,4阶格式：自行推导 （利用小程序）,非均匀网格：,问题： 会产生大量非物理解 （例： 1000个网格点算出1000个特征值） ； 可通过不同网格的对比，进行筛选,Copyright by Li Xinliang,17,离散化后得：,Step 2: 求广义解特征值问题（1） 即w为何值时（1）有非零解,(1),方法1： 全局法 一次计算出全部特征值,常用Q-Z分解法; 其他： 幂法、反幂法、Jacobi法，Householder ),狭义特征值问题,广义特征值问题,求解特征值是计算数学的主要研究方向，有大量成熟的方法,可借助软件包或Lapack库等 （自行到网上搜索）,通常关心最不稳定的扰动波,最大的那个,Copyright by Li Xinliang,18,方法2： 局部法,令：,方法： Newton法， 弦位法，抛物线（Muler）法,Newton 法：,弦位法：,差分化,抛物线法：,已知：,可连成一条抛物线,令：,求出新的值,消元法计算行列式（利用5对角特征）,计算出 后，求解方程（1）就可得到特征向量。,（1）,方程有奇性，可补充一个条件，例如给定某个点的值,Copyright by Li Xinliang,19,效果更好的方法 Malik提出的紧致差分格式求解 参考文献: 周恒等： 流动稳定性 （p. 10-13），国防工业出版社,非均匀网格的超紧致格式（4阶精度）,优点： 紧致网格基仅2点 精度高 4阶 直接适用于非均匀网格无需坐标变换,原方程组：,形式变换 变为一阶方程组,通常的做法,推导仓促 可能有误 请仔细推导,与y无关,Copyright by Li Xinliang,20,令：,带入差分格式：,得：,Copyright by Li Xinliang,21,令：,边界处表达式 (C1,D1,CN,DN) 根据边界条件而定,内点的表达式,特点： 离散形成的代数方程组呈 块两对角特征,求行列式，特征值都非常便利！,其余步骤与前文相同 可用矩阵广义特征值理论计算全部模态 也可利用 计算单个模态,计算域可以取为-1,1。 也可取为-1,0， 在中心线给对称（或反对称）边界条件,Copyright by Li Xinliang,22,边界条件,1） 如果取完整计算域 -1,1,2) 如果取一半计算域-1,0,处可设定对称或反对称边界条件,如设定对称条件，只能计算出对称扰动模态 如设定反对称条件，只能计算出反对称模态,对于槽道流，最不稳定模态是对称模态,但有些情况下，稳定模态在转捩过程中也发挥作用（例如感受性过程，见Zhong et al. JFM 556,55-103,2006）,按照内点方法计算,先根据处理内点的方法，求出所有点上的C,D值，再对边界点进行特殊处理（D的前两行设为0， C的前两行设为(1,0,0,0)及（0，1，0，0）,Copyright by Li Xinliang,23,具体解法（局部法）,求解,显然,因此只需计算每个4*4矩阵的行列式即可 （可直接写出表达式，也可用消元法计算）,编制好计算 的子程序后，可利用Newton法，弦位法，抛物线法等求解 ，得到复增长率,Copyright by Li Xinliang,24,计算出 后，利用消元法求解方程 ，即可得到特征向量,消元过程中，充分利用两对角块矩阵的性质，可减少计算量,独立消元,注： 由于方程有奇性，消元后，最后一个方程为0=0, 舍弃最后一个方程，并令最后一个未知数为1 （ ），即可解出 显然 （a为任意常数）也是原方程的解，因此 可用某一值（例如 ）归一化。,最终，得到 条件下，波数为 的最不稳定扰动波：,扰动波型函数的分布,Copyright by Li Xinliang,25, 13.2 边界层转捩的预测方法,1. 经验公式法,转捩位置,Mach 6 钝锥边界层表面的摩擦系数分布 (Li et al. Phys. Fluid. 22， 025105， 2010. Li et al. AIAA J. 46（11），2899-2913，2008),x,摩阻或热流,转捩起始点（transition onset),转捩峰（transition peak),充分发展湍流区,球锥的转捩Reynolds数,边界层外缘的Mach数,动量厚度定义的转捩Reynolds数,Copyright by Li Xinliang,26,2. eN 方法,LST理论：,积分起始点,扰动波进入中性曲线后，开始增长，局部增长率为,eN理论：扰动波增长到eN倍，即发生转捩,N值需要由实验（或经验）确定，通常为811, 即扰动增长4个量级（10000倍）左右。 在不可压缩及航空领域（亚、跨、超）较为成熟。 在航天领域（高超声速），还有待检验。,不足之处： 未考虑扰动波进入中性曲线前的衰减过程，没考虑感受性过程。,他人的改进： 苏彩虹，周恒等 考虑衰减过程 C. H. Su, and H. Zhou, Science in China G, 52 (1):115-123 (2009).,Copyright by Li Xinliang,27,3. PSE （抛物化扰动方程）法,优点： 1） 无需平行流假设； 2） 可处理非线性(N-PSE),Step 1: 得到扰动的控制方程,已知解,线性化,L-PSE,N-PSE,Step 2: 假设扰动具有波动形式,振幅，沿x方向是个缓变量 （相对y方向而言）,Step 3: 带入扰动方程，得到振幅 的控制方程,“缓变量”是个很有用的概念，可用来简化方程,Plantdl边界层理论就是利用“缓变量”的概念进行简化的。,LST的 方程是一维的 PSE的 方程是二维的,Step 4: 利用缓变量性质，舍弃方程中的椭圆项（为高阶小量），得到抛物化的扰动方程,沿x方向推进求解 （类似时间方向的处理），计算量相当于一维问题。 （“抛物化”的优势）,非线性项的处理方法与谱方法相似,Copyright by Li Xinliang,28,4. 转捩模型法 （间歇因子模型）,实际粘性系数,层流粘性系数,湍流粘性系数 （由湍流模型给定）,湍流间歇因子 （0表示纯层流，1表示纯湍流）,方法 1） 根据经验公式，给定 沿流向的分布 方法2） 给出 的发展方程，进行求解,Copyright by Li Xinliang,29,常识： 湍流的间歇性,外间歇性： 层流及湍流交替出现的现象,Mach 6 钝锥边界层的密度分布 （Li et al. PoF 22, 2010),边界层有清晰锐利的界面（层流区、湍流区“泾渭分明”）,湍流信号,层流信号,内间歇性： 湍流脉动的概率密度分布偏离Gauss分布（随机分布）,概率论（中心极限定理）独立