基于喷流拟序结构预测的SGS 模型比较研究1)
2021-11-09刘琪麟赖焕新
刘琪麟 赖焕新
(华东理工大学机械与动力工程学院,上海 200237)
引言
拟序结构普遍存在于湍流中[1-3],并影响湍流的发生、输运与耗散等过程,对湍流拟序结构的研究、对于理解湍流机理、分析流动相互作用有重要的意义.从大涡模拟来看,正确预测湍流拟序结构的关键就在于算法和亚网格尺度(sub-grid-scale,SGS)模型的准确性.
在算法方面,传统的大涡模拟方法求解网格滤波的控制方程,并采用SGS 模型模化亚网格尺度流动与网格尺度流动之间的相互作用.虽然近年来部分学者利用数值离散的特点,也发展出不需要添加模型的隐式大涡模拟方法[4-5],但为了针对性地讨论应用中亚网格模型的选择问题,本文仍采用高精度有限差分结合亚网格模型的大涡模拟经典算法,并讨论SGS 模型的性能.对大涡模拟算法更全面的介绍可以参考Sagaut 等[6]的综述.
目前文献中广泛出现的亚网格尺度模型是经典的常系数Smagorinsky 模型(Smagorinsky model,SM),已有众多文献报道该模型过度耗散了脉动动能[7-9].另一类是基于Germano 恒等式和Lilly 最小二乘法来动态调整模型系数的动态Smagorinsky 模型(dynamic Smagorinsky model,DSM),这类模型的局限性在于理论上要求在流场均匀的方向上计算平均的模型系数,而这一平均过程缺少物理上合理的解释,且造成分布式并行计算中不同处理器之间数据的频繁调用与交换,降低计算效率.从计算的角度,较为理想的模型是局部模型,计算效率高且能较好地克服SM 对湍流黏性估计过高的缺点.Piomelli和Liu[9]提出了局部动态模型(localized dynamic Smagorinsky model,LDSM),去除了DSM 对流动中需存在各向同性方向的要求,而是利用流动的历史信息来避免平均过程,计算效率得到提高.Lenormand等[10]提出的选择多尺度模型(selective mixed-scale model,SMSM),在边界层流动[10]、翼型绕流[11]和空腔流动[12]的预测中具有良好应用.这一模型借助涡矢量方向的变化来判断流动的状态是否发生了转捩,同时反映湍流的间歇性.Kobayashi 等[13-14]则利用流场中黏性耗散与拟序结构的关联性,提出了拟序结构模型(coherent-structure Smagorinsky model,CSM),并应用于各向同性湍流、湍流通道流以及横向射流的计算中,其模型简单且具有类似于动态模型的性能.Koyabashi[13]同时改进动能模型,提出拟序结构动能模型(coherent-structure kinetic-energy model,CKM).这些局部模型都有望克服传统模型的不足,但对于预测可压缩喷流中的拟序结构,目前尚未看到有关它们预测性能的分析与评估的文献.
典型的喷流拟序结构有流向涡、高低速流体间的条带[15]、剪切层边沿的卷吸流动结构[16]等,表现为长时间保持一定形状的大尺度质量团,且在时间和空间上具有强相关性[17].Lumley[18]最早提出采用本征正交分解(proper orthogonal decomposition,POD)方法提取拟序结构,且在自由剪切流、分离流、边界层转捩以及湍流通道流中得到广泛的应用.POD 方法的目标是针对物理场寻找一组空间正交的基函数,并按照基函数对物理量方差的贡献排序,例如对于脉动速度场的研究就是基函数对湍动能的贡献[19-20].本征正交分解的特点使得只需要少量的POD 模态就能够提取出流场中能量占主导地位的流动模式.此外,POD 分解的基函数直接由湍流场确定,不同于傅里叶变换或者小波变换需要人为指定基函数,因此表征的拟序结构更为客观[17].在POD 算法方面,虽然快照POD 算法能够显著减少计算量,但是有部分高阶模态没有参与计算.基于矩阵奇异值分解(singular value decomposition,SVD)的经典POD 算法包含完整的模态信息,比快照POD算法对舍入误差的鲁棒性高[19-20].此外数学软件Matlab 中提供了SVD 的库函数,也便于编制后处理程序.
本文对亚声速可压缩平行喷流进行大涡模拟.在验证计算结果的基础上,研究SM,CKM,SMSM,LDSM 和CSM 模型的特性,分析各个模型预测的瞬时涡结构.同时采用基于SVD 的POD 分解方法提取喷流拟序结构并分析各模型预测结果的差异,为可压缩喷流的大涡模拟提供参考.
1 数值计算与方法
1.1 亚网格尺度模型
控制方程组为Favre 质量加权滤波的三维非稳态可压缩N-S 方程组[21].方程中的亚网格尺度未知量由亚网格模型给出,包括亚网格尺度应力τij
以及亚网格尺度热通量
其中ρ,T,νsgs以及ui(i=1,2,3) 分别为密度、温度、模化的亚网格黏性系数以及速度矢量的3 个笛卡尔分量.γ 是量热完全气体的比热比,空气取γ=1.4.M为马赫数,湍流普朗特数常数Prt=1[22],而,顶标“-”表示网格滤波,“~”表示Favre 质量加权的网格滤波,所有计算变量均为网格滤波后的变量,因此除非特殊说明以下省去顶部符号.
Smagorinsky 模型(SM)中的亚网格尺度应力为
δij=0.对于喷流,取模型系数C1=0.17[23].τkk的模化根据文献[7]的建议得出.其中=Li j−δi jLkk/3,顶标“^”表示测试滤波运算,Ct为上一时刻的模型系数,初始值可以根据SM 模型取Ct=C1,当前时刻的模型系数C2是未知量,且式(5)对C2是超定的.C2的估计值由式(5)两边取αij的缩并得出[9],如式(6)所示.LDSM 模型的亚网格尺度应力与亚网格黏性系数通过将式(4)中的替换为得出
选择多尺度模型(slective mixed-scale model,SMSM)是根据局部涡矢量方向的脉动来调整亚网格黏性系数
其中C3=0.06,α 取0.5[24].为了反映湍流的多尺度特性,引入小尺度脉动的动能=(1/2),其中脉动速度=−.进一步,对网格滤波尺度的瞬时涡量进行测试滤波,得出.和两个向量之间的夹角 θ 反映了涡运动的脉动强弱.对于均匀各向同性湍流,在夹角 θ 的概率密度函数的峰值位置处有θ0=20°[24],据此,可判断局部流动的状态.模型引入了选择函数fθ0,当 θ> θ0时,取fθ0=1,表示涡方向在两个滤波尺度上变化较大,是无规则的完全湍流状态;否则,当 θ<θ0时,取fθ0=tan(θ/2)/tan(θ0/2),意味着湍流逐步减弱,亚网格黏性也应相应减弱,直到在层流区仅保留分子黏性.因此,SMSM 模型的亚网格黏性系数为
拟序结构模型是利用拟序结构与湍流中耗散的相关性来调整亚网格黏性系数vsgs=CΔ2|S(u)|,模型系数C=C4|FCS|3/2中常数C4=0.05[13],FCS=Q/E为拟序结构函数,其中
为速度导数Dij=∂ui/∂xj的第二不变量,它主要由Di j中反映当地旋转的反对称部分Aij和反映当地变形的对称部分Bi j组成,而为Di j的范数.其中
基于直接数值模拟(direct numerical simulation,DNS),流场的相关性分析表明,Q的绝对值与拟序小尺度涡附近能量耗散的强弱存在关联[13],因此|FCS|反映了流场中黏性与耗散的局部变化,而|FCS|的幂次“3/2”是根据不可压缩壁面流动中C∝y3,Q∝y2以及E∝ 常数(y是壁面的法向)而确定的.因此,CSM 模型的亚网格黏性系数为
拟序结构动能模型是在拟序结构模型的基础上,在亚网格黏性系数中显式地包含亚网格尺度动能ksgs,并做量纲调整后得到
至此,可以通过式(2)、式(4)来计算亚网格尺度热通量和亚网格尺度应力.虽然文献[12]指出亚网格尺度应力中各向同性的部分 τkk相对于热力学压力可以忽略,但根据Vreman[7]的建议,在基于SM 的LDSM 以及CSM 模型中仍然模化 τkk.
1.2 计算对象与数值离散
本文采用Hu 等[22]的亚声速平行喷流为考核算例.喷流出口雷诺数ReJ=DJρ∞UJ/µ∞=2000,马赫数MJ=UJ/c∞=0.9,其中DJ=1、UJ=1 分别为归一化的喷口宽度和喷口速度,参考值ρ∞,µ∞和c∞分别为来流密度、动力黏度以及环境声速.计算域与网格划分如图1 所示,流向(x)、横向(y)与展向(z)的尺寸分别为Lx×Ly×Lz=24×15×3,采用笛卡尔网格进行离散,对应方向的结点数分别为 181×181×16.为了提高近喷口区域和喷流唇线附近的网格分辨率,采用双曲正切函数对x和y方向的网格进行拉伸,而喷流在z方向采用周期边界条件,划分均匀的网格.基于Celik 等[25]提出的大涡模拟分辨率指标,由于课题组前期工作已证明该网格能够分辨湍动喷流中75%以上的湍动能,考虑到喷流的雷诺数较低,该计算网格可以满足大涡模拟的网格质量要求[26].
图1 计算域与网格(3 个方向都每隔3 条网格线画一条线)Fig.1 Computional domain and the Grid(One in every three lines is shown in each direction)
采用有限差分方法开展数值模拟,本文的计算源程序以Sandham 课题组[27-28]的(shock-boundary layer interaction,SBLI)代码为基础,计算精度和可靠性经过了充分的考核验证[22,29].在此基础上,本文进一步采用OpenMP 并行化,使该代码可用于内存共享式的小型计算工作站.
空间离散的1 阶与2 阶导数均采用5 点4 阶中心差分格式计算,边界点采用Carpenter 提出的单侧4 阶稳定差分格式[30],从而保持空间离散的整体高精度.为了抑制差分计算引起的非物理振荡,算法上采用了熵分裂方法,通过将对流项分裂为守恒形式与非守恒形式的加权求和来实现[31].同时,在x和y方向的计算边界上,采用无反射边界条件来抑制非物理振荡波在计算域边界上的反射,算法上通过保留走出计算域的特征波而将走入计算域的特征波幅值置零来实现[32].在z方向上采用周期边界条件来模拟平行喷流在展向的无限延伸.
时间离散采用低存储显式3 阶龙格库塔算法,这一算法完成3 个子时间步的推进仅需要两个存储流场的数组,对内存资源的消耗低[21].为便于模型特性的研究比较,时间步长取定值Δt=0.003 075,且满足数值稳定性条件CFL<1.在初始条件方面,给定喷口平面上流向速度u、温度T以及压强p的形线,并用喷口平面上的物理量来初始化全场.喷口的速度分布采用双曲正切函数给定,并施加伴流速度U0=0.1,可得
在数值结果方面,虽然DNS 能够准确给出流场结构的信息,但是DNS 对计算资源的消耗量高,而且流场后处理的数据量大.以文献[22-33]中提供的DNS 数据为基础,课题组前期工作已开展了计算代码的考核[34-35]、LES 的网格分辨率分析[26]以及局部模型的研究[21,36],算法的可靠性与模型都已得到验证.本文在已有工作的基础上,进一步分析局部模型揭示的喷流拟序结构.
1.3 POD 方法
POD 方法采用一组正交基函数 ψ(x) 的线性组合来表示物理场 φ(x,t),组合系数为c(t),即
所得的基函数也称为POD 模态.POD 分解的目的是使这种表示的误差在平方平均的意义下最小,并且在给定误差下最小化表示所需要的模态个数,从而得到最佳的基函数[20].这组最佳基函数能够最优地表征原物理量的方差,对于脉动速度场,这种最优是模态反映的湍动能随模态阶数增大而迅速的减小[19],脉动速度的第一阶POD 模态对湍动能的贡献最大.详细的理论介绍可参考文献[20].
在计算中,可对流场的速度分量进行时间采样,得到关于标量 φ 的N个流场样本,对于网格点数为M=Mx×My的二维样本,可将每个样本都写成M×1的列向量,由此得到M×N的样本矩阵
在离散求解时,寻找最佳基函数的过程在数学上等价于求样本矩阵F的SVD 问题[20].采用Matlab 中低存储的“svd()”库函数计算F的奇异值分解
其中,左矩阵 Ψ 的各列是正交的,Ψ=(ψ1,ψ2,···,ψN),它的第i列 ψi即为物理场的第i阶POD 模态,右矩阵C的各列也是正交的,C=(c1,c2,···,cN),它的第i列ci即为模态 ψi对应的系数,反映POD 模态与时间有关的特征[19],而S为奇异值矩阵,S=dig(s1,s2,···,sN),S的对角元为从大到小排列的奇异值si(i=1,2,···,N),其他矩阵元为零,λi=代表了与第i阶POD 模态对应的能量,且对各阶模态可以定义能量贡献率
以及累计能量贡献率
另外,由于各阶POD 模态 ψi是相互正交的,所以λ衰减越快,当前模态所含有的能量相比下一阶模态的差就越大,当前模态所表征的物理场的相关性就越强.
2 SGS 模型的特性
2.1 平均流场
在分析湍流统计特征时,为排除初场对统计结果的影响,取无量纲时间t=48~ 216 的流场数据,对应流体以喷口速度流过2~ 9 倍的计算域长度,以计算时间步长Δt为采样间隔,得到54 640 个流场样本.图2 以LDSM 模型为例,显示了在唇线上无量纲位置x=2,6,8,10,20 处的流向脉动速度的功率谱,谱线具有一段固定的斜率且能观察到惯性区,反映了湍流的能量级串现象,也说明流动已经达到了完全湍动状态.
图2 LDSM 模型预测的脉动速度无量纲功率谱Fig.2 Dimensionless power spectra of velocity signals predicted by the LDSM
喷流主流流体与环境流体存在速度差,形成剪切层.图3 所示为喷流剪切层的发展.考虑喷流中心线速度Uc与伴流速度U0,图3(a)采用速度半值宽δ0.5,即流向速度u= 0.5(Uc+U0) 时的横向位置y,显示出剪切层的发展呈现两个线性阶段,即初期的缓慢增长段与下游的快速增长段,这对应喷流势流核末端x≈6 处的流动转捩现象.对比5 种模型的结果,SM 模型预测的剪切层横向发展程度较其他模型大,表明该模型过度耗散的特点.图3(b)以LDSM 模型为例,在相似坐标中绘出了喷流不同流向位置的速度形线,相似坐标的横向位置采用速度半值宽归一化,平均喷流速度差〈u〉−Uc采用主、伴流的速度差ΔUc=Uc−U0归一化,〈·〉 表示时间平均.图3(b)中,速度形线在入口呈平顶状分布,而在下游则快速发展为自相似的速度形线,可以看到x=4 时,自相似形线还没有建立,但当x≥8 时,速度形线与实验喷流[37-39]的自相似数据点吻合良好,也验证了计算结果.
图3 喷流剪切层的发展Fig.3 Developments of the jet shear layer
喷流流体与环境流体的强烈掺混伴随着能量的耗散.对比5 种SGS 模型的亚网格尺度平均黏性耗散率,定义为.图4 为在展向中心平面(z=1.5)上的分布.由于喷流入口施加的平均来流条件,流体喷出后在近喷口区基本为层流状态,随后在势流核末端(x≈6)附近发生流动转捩,并在下游发展为充分湍流.图2 中,近喷口的脉动水平比势流核及其下游的区域低1~ 2 个数量级,进一步表明来流的层流特征,此处与湍流脉动相关的亚网格尺度黏性耗散应趋近于零.但在图4(a)中,SM 模型预测的在近喷口区内呈两条平行的条带,表明SM 模型对层流扰动敏感,不能正确反映湍流黏性耗散率的分布区域.图4(b)~ 图4(e)分别为CKM,SMSM,LDSM 以及CSM 的预测结果,这4 种局部模型虽然能够正确地反映近喷口区的分布特征,在转捩区(x=6~9)内呈现一对峰值,但是以CSM 模型为参考,CKM,LDSM,SMSM 以及SM 模型预测的峰值耗散不同,分别为CSM 模型预测值的3.0,1.8,1.3 和0.9 倍.其中CKM 预测的峰值接近其他模型的1~ 2 倍,稍后分析的分布对转捩区内涡结构的影响.
图4 对比在中心平面(z=1.5)上的分布Fig.4 Comparisons of the distributi on of ?in the center plane of the jet (z=1.5)
2.2 瞬时流场的涡结构
分析无量纲时间t=207 的瞬时流场,此时喷流以喷口速度流过8 倍以上的计算域长度,已经充分排除初场的影响.参考Hu 等[22]的部分DNS 结果,以喷流中心线速度〈u〉 衰减为0.99 时的流向位置得出平均的势流核长度为x=6.43.图5 为Q=0.1 的等值面,图5 中红色切片位于Hu 等DNS 结果[22]的势流核末端.在预测流动转捩方面,图5(c)~ 图5(e)中,SMSM,LDSM 以及CSM 模型都能预测大尺度展向涡在平均势流核末端附近的破碎,以及下游区域内湍流多尺度涡结构的生成.然而对于SM 模型,图5(a)在切片位置下游仍然能够观察到流向的大尺度涡结构(如图5 中箭头所示),结合图4(a)可知,湍流中小尺度涡的发展可能受到近喷口过度的亚网格尺度耗散的抑制,这与文献[7-9]对SM 模型过度耗散特点的报道一致.图5(b)中,CKM 模型的预测结果与SM 类似,同样可能由于该模型中局部亚网格尺度耗散较强.图6 为展向涡 ωz的等值面,同样采用红色切片标记DNS 结果[22]的势流核末端位置.可以看到,ωz的等值面在势流段内呈平行的涡片,随着流动向下游的发展,ωz等值面在势流核末端附近快速卷起,随后涡结构被拉伸并在下游发展出复杂的三维结构.图6(a)中,SM 模型预测的涡卷起位置位于切片之后,涡片偏长,对应图4(a)中近喷口区的不合理耗散,SM 模型抑制了流场中展向涡的发展.图6(b),图6(d),图6(e) 中,CKM,LDSM 以及CSM 模型在势流核末端位置附近都能够观察到展向涡的卷起(如图6 中箭头标记)以及下游的涡破碎,但LDSM 与CSM 模型在切片后预测的小尺度涡结构比CKM 模型丰富,再次表明CKM 模型的局部强耗散可对小尺度涡结构产生抑制.图6(c) 中,SMSM 预测的涡卷起不显著,但切片下游的预测结果与LDSM,CSM类似.相比SM 与CKM 模型,SMSM,LDSM 以及CSM 模型更能反映势流核下游展向涡的多尺度特性.
图5 对比 Q=0.1 的等值面,采用流向速度着色,红色切片位置标记了文献[22]中DNS 结果的势流核末端:x=6.43Fig.5 Comparisons of the iso-surface of Q=0.1,colored by streamwize velocity.The red slice marks the end of the potential core x=6.43 of the DNS results in Ref.[22]
图6 对比 ωz 的等值面,红色切片位置标记了文献[22]中DNS 结果的势流核末端:x=6.43Fig.6 Comparisons of the iso-surface of ωz,the red slice marks the end of the potential core x=6.43 of the DNS results in Ref.[22]
2.3 速度场的POD 模态分析
对5 种SGS 模型的预测结果进行POD 分解.平行喷流主要表现为x-y平面内的二维流动特征,因此在湍流统计的流场数据中,截取三维流场的中心平面(z=1.5)为样本,以10Δt(=0.030 75) 为采样间隔,得到4685 个样本.针对5 种SGS 模型,提取脉动速度分量u′(=u−〈u〉),v′(=v−〈v〉) 和w′(=w−〈w〉)的POD 模态.脉动速度u′,v′和w′的第一阶模态分别对原始脉动速度场的脉动强度:〈u′u′〉,〈v′v′〉 和〈w′w′〉的贡献最大,并在这一意义下脉动速度的第一阶模态与原始脉动速度场最为接近.图7 为流向脉动速度u′模态的累计能量贡献率 η,可以看到,累计能量贡献率曲线在前200 阶内就达到了90%以上,这200 阶模态仅为总模态数的4%,因此曲线 η 反映了POD 的快速收敛特性.但5 种SGS 模型的 η 具有不同的增长速度,SM 模型最快,而SMSM 模型最慢.脉动速度v′与w′的 η 曲线具有与图7 相似的快速收敛特征,因此不再绘出.图8 为脉动速度分量在前20 阶模态的能量贡献率衰减曲线 ξ.可以看到,5 种SGS 模型得出的 ξ 曲线都在前5 阶内迅速衰减,说明低阶模态占主导且模态表征的流场具有强的相关性.对比图8 中第1 阶模态的能量贡献率,图8(a)中u′模态的衰减曲线显示SM 模型的第1 阶模态贡献率最高,SMSM 模型最低;图8(b)中横向速度脉动v′模态的衰减曲线显示CKM 模型的贡献率最高,而SMSM 模型依然排最低;图8(c)展向速度脉动w′模态的衰减曲线中SMSM 模型的能量贡献率排序有所提高.图8 表明,第1 阶模态的能量贡献率已经体现出了不同模型之间的差异,接下来针对第1 阶模态展开分析.
图7 u′ 模态的累计能量贡献率的收敛曲线Fig.7 Convergency curves of the cumulative energy contribution for the modes of u′
喷流主要表现为流向的动量传递.对比5 种SGS 模型预测的流向速度脉动u′的第一阶POD 模态,图9 显示了模态幅值在中心平面(z=1.5)上的分布,采用接近零的模态幅值等值线u′=±0.003 来显示模态的轮廓.图9(a)中,SM 得到两条沿流向的带状模态,且在喷流下游x=12~24 的范围内一直保持高的模态幅值,对应图5(a)中大尺度流向涡结构对流场的搅动.图9(b)中,CKM 模型预测的模态也在下游保持沿流向的高幅值带状模态,在转捩区(x=6~9)内的模态强度呈现一对低谷区,位置恰好匹配图4(b)中的峰值耗散区位置.这表明流向速度脉动的POD 模态对亚网格尺度的黏性耗散率敏感,局部的强耗散导致脉动强度〈u′u′〉 的降低,从而改变了POD 模态的形状.图9(d)、图9(e)分别为LDSM与CSM 模型预测的模态,两者都沿流向两侧发散开来,并且在喷流下游呈现分岔或破碎的形态,模态幅值也逐渐减小,反映下游的多尺度结构以及脉动能量的逐渐消散.图9(c)为SMSM 模型预测的模态,模态轮廓在流向上呈断续的条状,局部的模态幅值低,这可能与图8(a)中SMSM 模型第一阶POD 模态的能量贡献率低有关.而在下游,SMSM 的模态呈现与CSM,LDSM 模型预测结果相似的分岔形态.因此,与SM 和CKM 模型相比,CSM,LDSM 以及SMSM 模型都合理地反映了u′模态的形态.
图8 脉动速度分量的POD 模态的能量贡献率Fig.8 Energy contribution rate of POD modes for the fluctuate velocity components
图9 流向脉动速度 u′ 的第一阶POD 模态云图,实、虚等值线分别表示 u′=+0.003和u′=−0.003Fig.9 Contours of the first-order POD modes of u′,with solid contour lines for u′=+0.003 and dashed contour lines foru′=−0.003
图10 为横向脉动速度v′的第一阶POD 模态,等值线v′=±0.003 显示出模态的轮廓.同时,图10中采用面内矢量 (u′,v′) 的POD 模态来显示喷流中的流动模式.观察图10,v′的POD 模态呈正负相间的肋状,并沿流向排列.在肋状模态上下端部的附近,矢量箭头显示出环形的流动模式,在喷流上部的红圆标出了逆时针环流,下部的蓝圆标出了顺时针环流.环流模式穿过了主流与环境流体,在圆环的下游侧,高速的主流流体被送入两侧的低速环境流体中,而在环流的上游侧,低速的环境流体被带入主流流体中.同时,v′模态轮廓的横向尺寸也沿着流向逐渐增大.可见,(u′,v′) 的POD 模态反映了喷流的流动卷吸现象[16].对比图10(a)~ 图10(e),SM,SMSM,LDSM 以及CSM模型预测的模态都反映了喷流的卷吸现象,其中SMSM 模型预测的旋涡尺度小,这与图8(a) 和图8(b) 中的低模态贡献率对应.而图10(b)中,CKM模型预测的模态基本没有反映该区域内的流动卷吸,结合图4(b),脉动速度分量u′和v′可能受到CKM模型在x=6 附近局部的强耗散抑制,导致未能预测出明显的流动卷吸.
图10 横向脉动速度v'的第一阶POD 模态云图与面内矢量 (u′,v′) 的第一阶POD 模态矢量图(实、虚等值线分别表示v′=+0.003和 v′=−0.003,红圈标示逆时针环流、蓝圈标示顺时针环流)Fig.10 Contours of the first-order POD modes of v' and the in plane vector (u′,v′)(Solid contour lines for v′=+0.003 and dashed contour lines for v′=−0.003.Red circle marks counterclockwise circular flow and blue circle indicates clockwise circular flow)
图11 为展向脉动速度w′的POD 模态云图,等值线w′=±0.003 显示出模态的轮廓.图11(b)~图11(e)分别为CKM,SMSM,LDSM 以及CSM 模型预测的结果,w′的POD 模态在转捩区(x=3~9)内呈现规则的脊状排列,红色表示脉动速度w′指向平面外,蓝色表示脉动速度w′指向平面内,强弱相间的分布反映了w′模态对涡结构的展向拉伸,如图11(f)所示.图11 中的脊状模态显示涡拉伸主要发生在转捩区,而在下游的模态轮廓呈碎块状,强度也逐渐减弱,对应涡结构逐渐发展为三维结构,展向的拉伸现象不再显著.图11(a)中,SM 模型预测的模态没有反映转捩区内的脊状排列,在下游呈大团块,位置与图4(a) 中SM 模型喷口附近的过度耗散区域相对应.图11(b)中,CKM 模型预测的模态虽然呈现了脊状排列,但是模态轮廓在下游x=21 处仍然呈大尺度块状,与SM 模型的预测结果相似.相比而言,SMSM,LDSM 以及CSM 模型的预测结果则能够合理地反映展向拉伸的流动模式.
图11 展向脉动速度 w′ 的第一阶POD 模态云图Fig.11 Contours of the first-order POD modes of w′,with solid contour lines for w′=+0.003and dashed contour lines forw′=−0.003
图11 展向脉动速度 w′ 的第一阶POD 模态云图 (续)Fig.11 Contours of the first-order POD modes of w′,with solid contour lines for w′=+0.003and dashed contour lines for w′=−0.003 (continued)
在计算效率方面,表1 对比了5 种SGS 模型的CPU 时间,以SM 模型为基准,LDSM,SMSM,CKM,CSM 的计算时间分别为SM 模型的:1.64,1.28,1.22 和1.00 倍,其中LDSM,SMSM 以及CKM 模型需要进行滤波运算,因此需要更多的CPU 时间,CSM 模型的计算时间最接近SM 模型.
表1 每时间步内每网格点所占用的CPU 时间Table 1 CPU time per grid point per time step
3 结论
本文开展了亚声速可压缩平行喷流的大涡模拟.在分析平均流动与耗散的基础上,针对瞬时涡流以及脉动速度POD 模态表征的拟序结构,对SM,CKM,SMSM,LDSM 和CSM,5 种SGS 模型的性能进行研究,得出以下主要结论:
(1) 初步揭示了u′模态从流向的带状到下游分岔或破碎形态反映的多尺度特征、v′模态沿流向的肋状排列与尺寸增长、矢量 (u′,v′) 模态的环流模式表征的流动卷吸、以及w′模态表征的流场结构在展向上受拉伸的模式.
(2) SMSM,LDSM 以及CSM 均较好地反映出湍流的多尺度特性,清晰地分辨了小尺度结构的发展过程.而SM 与CKM 模型则未能有效地预测湍流中的小尺度涡结构及部分流场特征模态.
(3) 脉动速度场的POD 模态对SGS 模型的亚网格尺度耗散敏感,表现为CKM 预测的峰值耗散区对应u′模态的低谷区因而环流模式不显著,以及SM 未能反映转捩区内w′模态的脊状拉伸模式.SMSM,LDSM 和CSM 克服了以上模型的不足之处,其中CSM 模型同时兼有计算效率较高的优势.