用PIV数据估算槽道内湍流动

范文一：用PIV数据估算槽道内湍流动能耗散率

　第55卷　第7期　　　　　　　　　　　　化　　　工　　　学　　　报　　　　　　　　　　　　　Vol.55　№7　　2004年7月　　　　　　Journal　of　Chemical　Industry　and　Engineering　(China)　　　　　　July　2004

研究论文

用PIV数据估算槽道内湍流动能耗散率

王汉封　柳朝晖　郭福水　郑楚光

(华中科技大学煤燃烧国家重点实验室,湖北武汉430074)

摘　要　湍能耗散率的准确获取对工程实际问题和湍流的理论研究都有着重要的意义.但由于湍能耗散率定义

的固有复杂性,其测量一直都是一项有相当难度并具有挑战性的工作.运用PIV对水平槽道内湍流流动进行了测量,详细介绍了运用大涡模拟中亚网格(SGS)应力模型估算湍流耗散率的大涡PIV方法,分析了运用该方法估算湍流耗散率的优势.将各种不同方法所得到的结果与DNS结果进行了比较,发现运用修正的Smagorinsky模型所得结果与其符合最好.

关键词　湍能耗散率　PIV　大涡模拟　亚网格模型　湍流测量中图分类号　TQ021.1　O357.5　　　　文献标识码　A

文章编号　0438-1157(2004)07-1066-06

ESTIMATIONOFTURBULENTKINETICENERGY

DISSIPATIONRATEINCHANNELFLOWBYPIV

WANGHanfeng,LIUZhaohui,GUOFushuiandZHENGChuguang

(NationalLaboratoryofCoalCombustion,HuazhongUniversityofScienceand

Technology,Wuhan430074,Hubei,China)

Abstract　Accurateestimationofturbulencedissipationrateisveryimportantforengineeringandacademic

studyofturbulence.Itischallengingtomeasurethedissipationrateduetoitscomplexdefinition.ThemeasurementoftheturbulenceflowinahorizontalchannelwasperformedbyusingPIV,andthusthemethodofusingSGSmodeltoestimatetheturbulencedissipationratewasintroducedandtheadvantageofsuchmethodwasanalyzed.AlltheestimationresultswerecomparedwiththeDNS'sresults.TheoutcomeofthecorrectedSmagorinskymodelfittedtheDNS'sresultsbest.

Keywords　turbulencedissipationrate,PIV,LES,sub-gridscalemodel(SGS),turbulencemeasurement

量纲分析的方法进行封闭,是各种湍流模型中最不

精确的部分.另一方面,由于湍流的多尺度特性,常被应用于化工、燃烧等工业过程中,而湍流耗散率又对这些过程有着重要的影响.因此,准确的测量湍能耗散率不仅可以对评价及改进湍流模型提供依据,而且对许多工业过程的设计也有一定的意义.

引　言

湍流条件下流动的能量方程中会出现湍能耗散项,该项代表流动中湍动能到热能的转化速率.由于耗散项定义的复杂性,对湍能耗散率的研究一直都是湍流理论研究中一个较为薄弱的环节.在湍流模型研究中,耗散率输运方程中的耗散项大多采用

　　2003-05-12收到初稿,2003-09-05收到修改稿.

联系人:柳朝晖.第一作者:王汉封,男,28岁,博士.基金项目:国家自然科学基金项目(No.50276021)及国家重大基础研究规划项目(No.G199902207)共同资助.

　　Receiveddate:2003-05-12.

Correspondingauthor:Dr.LIUZhaohui.E-mail:zliu@mail.hust.edu.cn

Foundationitem:supportedbytheNationalNaturalScienceFoun-dationofChina(No.50276021).

　第55卷第7期　　　　　　　　　王汉封等:用PIV数据估算槽道内湍流动能耗散率

·1067·

湍流耗散率的传统测量方法都是由热线风速仪(HWA)的测量数据处理得到.运用HWA测量湍流耗散率可以依据其实现的方法不同分为三类:从时间序列获得能谱、空间关联系数、直接测量脉动速度的空间梯度等不同方法[1].在能谱法中为了得到较好的一维能谱E1(k)的解析值,必须拥有很高的速度采样频率,该方法尤其适合用于采样频率较高的热线风速仪,通常采用的一维能谱的关系式中已包含了Taylor假设.Azad和Kassab[2]采用该方法得到边界层流动和充分发展管流内的耗散率分布情况.Azad[2]和Elsner[1]对热线长度l与Kol-mogorov尺度η的相对大小对该方法测量结果的影响作了研究,发现在l/η=5时,由于热线长度所带来的耗散率测量误差将在15%左右.所谓的空间关联系数法本质上与能谱法类似.该方法采用两个独立热线探针,通过改变两热线间距离,测量不同距离上的脉动关联系数,以重建得到湍流能谱[1].由于该方法实现较为复杂,文献中一般没有采用.速度梯度方法是运用两个或更多的探针同时测量,直接得到所需要的脉动速度的空间梯度,通过耗散率定义式(2),或简化式(3)得到湍流耗散率.Tsinober[3]运用专门设计的具有21根热线组合探针,通过式(2)得到了耗散率.An-dreopoulos与Honkan[4]采用了9根热线的组合探针,通过简化假设得到了充分发展管流内耗散率分布.这种方法虽然可以对耗散率定义式中各项进行直接测量,但是由于多热线探针自身的复杂性,使探针空间分辨率受到了限制.同时随着探针结构的复杂化,其本身对流场的干扰以及各热丝间的干扰也将给测量带来更大的系统误差.

PIV测量技术是最近10年来随着激光技术、图像采集技术以及计算机技术的发展而发展起来的一种新型的流动测量技术.相对于以往测量方法其最大的特点是实现了同时多点测量的“平面测量”.运用PIV技术,不仅能用已有的脉动速度空间梯度的方法进行湍能耗散率的测量,而且有更适合PIV测量特性的“大涡PIV”方法可以得到湍流耗散率[5].

本文将对上述各种方法进行简要介绍与比较,并详细给出大涡PIV方法的原理及实现步骤.采用脉动速度梯度法和大涡PIV等方法对矩形水平槽道内湍流耗散率进行估算,给出不同方法所得结1　耗散率定义

由于黏性作用,真实流体的流动总是具有耗散性的,这种由黏性引起的机械能到内能的转化称为机械能耗散率.机械能耗散主要集中在最小湍流涡团的尺度上[1],因此,在本文中将只研究湍能耗散率.

某一点的湍能耗散率能够通过其定义[3]得到,一般情况下,其表达式如下所示

uε+ xi xj xi

(1)

在三维情况下,式(1)可展开为

ε2

++++ x1 x1 x1 x2

u u u u u2++++2+

22333

2 x x+2 x x+2 x x

213132

(2)

运用局部均匀各向同性的假设,式(2)可以简化为如下的形式

u u++12ε=6 x x221

(3)

另外在各向同性的湍流中,耗散率可以简化为用长度尺度和均方根脉动速度来表示,即ε=15ν(u′/

λ),式中λ为Taylor微尺度.利用湍流平衡假

设,可以下式表示

ε=C(u′/l)

(4)

式中　u′是均方根脉动速度,l是积分长度尺度,C为比例常数.

2　PIV方法

2.1　普通PIV方法

PIV系统能够同时测量一个平面上的速度分布,因此,速度的空间导数和其他空间关联量从理论上说都可以直接从PIV测量的速度场中得到.

然而,值得注意的是PIV的空间分辨率是有限的,可用下式表示

Lr=(LIC/LCCD)Lv

(5)

式中　Lr是PIV系统可分辨的最小长度尺度,LIC和LCCD分别为查问区间大小和CCD的分辨率(通常用像素数表示),Lv表示CCD所实际拍摄范围的尺寸.显然对一特定的PIV系统,要求系统分辨率提高则必须减小CCD拍摄范围.例如在本实验

中,LIC=64像素,LCCD=1280像素,依据平均流,

·1068·

化　　　工　　　学　　　报　　　　　　　　　　　　　　2004年7月　

80μm,为了能够使PIV的分辨率达到这样的要求,PIV实际拍摄范围Lv只能小于1.6mm,这样的范围对于实际应用显然太小.因此,为了发挥PIV测量特点又避免上述情况发生,需要一种新的适用于整个流场测量的耗散率处理方法.

作为互相关算法的结果,PIV得到的速度矢量是查问区间范围内速度的空间平均值.因此,应用PIV测量速度场及其他湍流量时,必须考虑PIV所测量的速度与真实速度间的差别.可以将查问区间内的速度对中点处速度u作二元Taylor展开得到下式,式中u(x,y)代表查问区间中心点的速度.

u(x+h,y+k)=u(x,y)+h+k×

1!u(x,y)h+k

2! x y

的性质,是通过“亚网格模型”进行求解的.基于上面的讨论,可以通过计算大涡模拟中亚网格尺度上的湍能传递量来估算湍能的耗散率.下面将给出

大涡模拟的基本思路及计算耗散率的具体方法.瞬时速度和压力可以表示为

ui= Ui+ ui,pi= Pi+ pi

(9)

式中　U i,Pi是滤波后可分辨的速度及压力,即PIV可以测量的尺度范围内的量;u pi是速度i, 和压力被滤除掉的部分.对瞬时不可压流体的N-S方程进行滤波可以得到LES方程

Ui Ui (Ui Uj) U i τij

=0,+=-2-(10) xi t xj xj xj xj

式中　τij=uiuj-U iU j是亚网格应力张量,须用亚网格模型进行模拟.将上式乘U i,可整理得到

(6)

u(x,y)+…

湍能方程

k k k

+ Uj=- P Uj+v-τUi-ij jj xj

U i Ui

+2τSijij jj

(11)

在查问区间上对h和k积分,然后取平均,可以得到空间平均速度值,即PIV所得到的结果,式中Δ表示查问区间的尺寸.

u=u+

Δ2++2

24 x y2

Δ4+…++19210 x3 x y10 y

(7)

式中　 Sij( U xi+ U xj),是应变率张j/ i/ 2

量.上式中的最后一项表示从大涡模拟可解析的尺度传递到亚网格尺度的湍能,从前面的讨论可以知道,这一项就等于黏性尺度上的湍能耗散率.因此,湍能耗散率能够通过下式计算得到

ε≈〈ε〈τSij〉ij SGS〉=-2

(12)

上式右边第二项在数量级上与Δu/24l2相当,l

是积分长度尺度.由此可见,当满足Δ l时,真实速度u可以用PIV测量得到的速度U 来近似.同样,下面的关系式也都会得到满足,括号表示系综平均

〈ui〉≈〈 ui〉,u′≈ ui-〈 ui〉,u′u= ui uj-〈 ui〉〈 uj〉(8)ii′j

3　实验简介

本实验是在一小型循环风洞上进行的.实验段截面为0.3m×0.03m,长度为6.0m的水平矩形槽道.在距离突扩稳定段0.8m处设置了玻璃的测

量段.实验系统布置及PIV系统分别如图1、图2所示.示踪粒子采用卫生香燃烧生成的烟雾,从扩容室引入风道内,经过轴流风机及风道弯管部分的作用,已与空气十分均匀地混合.示踪粒子的平均粒径在5μm以下,其流动跟随性和光学散射特性

2.2　大涡PIV法

在Reynolds数较高的流动中,可以认为湍能是在与积分尺度相当的大尺度涡团上产生的,并随

着湍流涡团的不断破碎向较小尺度的涡团传递,最终在黏性尺度范围内被耗散掉.在产生湍能的大尺度和耗散湍能的黏性尺度之间,可以认为湍能既不产生也不耗散,在这个范围内的涡团仅仅是将由较大涡团传递来的湍动能继续传递给较小尺度的涡团,而所传递的能量与最终耗散掉的能量在数量上是相等的.

大涡模拟(LES)中所谓的截断尺度就是处在上述只有湍能传递的尺度范围之内,大涡模拟用截断尺度对流动的N-S方程进行过滤,对大于该尺度的运动直接进行求解,而对小于该尺度的运动,,Fig.1　Experimentalwind-tunnelset-up

1—valve;2—axialblower;3—windtunnel;4—elbowbend;

5—stabilizationsection;6—constrictionsection;7—;8—sudden

　第55卷第7期　　　　　　　　　王汉封等:用PIV数据估算槽道内湍流动能耗散率

·1069·

都很好.主要实验参数如表1所示

Fig.2　PIVsystem

1—lasersource;2—lightsheet;3—testsection;4—powersupply;5—CCD;6—synchronizer;7—PC

Fig.3　VelocityandfluctuationvelocitycapturedbyPIV

insky模型和梯度模型时,出现的所有与w方向有关的梯度量,也必须用已测量得到的u、v两方向

的相应量作为近似.并采用量纲近似法、速度梯度

Table1　Mainexperimentalparameters

Pulse

interval/μs10Meanvelocity/m·s-14.3

Resolution

/μm·pixel-1

5.6Reb4300

Interrogation

area/pixel64×32Reτ213

Samplingrate/Hz3.75

Lens

focus/mm105Kolmogorovlengthscale/μm

法,以及大涡PIV等不同的方法,对PIV系统得到的数据进行处理.

大涡PIV方法中的亚网格应力模型分别采用了修正后的Smagorinsky模型和梯度模型.其中修正的Smagorinsky模型如下

τSij=-CSA2Δ2 S Sij ij=vt

Integral

lengthscale/mm

1.5

(14)

为保证实验中的分辨率满足使用大涡PIV方法的条件,本实验中测量耗散率时使用放大率1∶1,在分辨率为5.6μm/pixel的状态下,PIV系统

连续采样1000帧,对水平槽道上壁面以下6mm的范围测量.

式中　 S =S ijSij;CS是Smagorinsky常数,取0.16;Δ是亚网格特征尺度,本文中取PIV查问区间的尺度;A为壁面修正系数,A=[1-exp

(-Y+/25)]1/2.梯度模型由Clark等于1979年提出.

τij=

Δ( U i/ xk)( U j/ xk)(15)

4　结果与讨论

图3中给出了分辨率为31.7μm/pixel时,

PIV测量结果与Moser等用DNS计算Reτ为180的结果对比.本实验中,uτ由槽道平均速度估算得到,对应的Reτ为213.所有的量都用壁面滑移速度uτ来量纲1化.由图中可以看出,PIV测量得到的结果与DNS结果符合得非常好,只是uu+在壁面附近有一定差别.这是由于PIV受到查问区间的限制,无法在十分贴近壁面的速度梯度很大的范围上准确测量,且在近壁区域内测量中错点偏多造成的.由图中还可看出,实验中非主流方向的脉动要明显弱于主流方向,因此,假设两个非主流方向的脉动特性相同,即w′=v′.要使用式(4)量纲分析方法时,其中表示脉动强度的均方根脉动速度u′可通过湍能来计算得到.

22222

k=(u′+v′+w′)/2≈(u′+2v′)/2

式中所有的量都是过滤后的值,该模型适合于处理

PIV测量得到的数据.

图4中给出了运用不同方法所得到结果的比较情况.所有的结果都用壁面滑移速度和运动黏性系数ν对其进行量纲1化.从三种方法所得耗散率分布可以发现,耗散率在壁面附近均出现峰值,但相对于速度梯度法和大涡PIV方法而言,量纲近似方法所得耗散率的变化较平缓,耗散率最大值以及最大值所出现的位置也均与其他两种方法以及DNS结果存在较大的区别.这是因为在量纲分析方法中所用到的湍流积分长度尺度l被假设为一个固定值,而l在壁面附近和槽道中心将有很大的变化,这一假设将带来较大的误差.图4中所给出的速度梯度法的结果,是依据近似式(3)统计得到的,而大涡PIV方法中所采用的是带有壁面修正的Smagorinsky亚网格模型和梯度模型.三种模型,(13)

·1070·

化　　　工　　　学　　　报　　　　　　　　　　　　　　2004年7月　

峰值,然后随着距壁面的增加而迅速下降,并在很大范围内保持恒定.由于PIV的分辨率仍远大于Kolmogorov尺度(约10倍),脉动速度的空间导数将无法准确得到,所以,在这样的分辨率下采用脉动速度空间梯度法所得结果将会有较大的误差

5　结　论

PIV技术能够同时测量一个平面上的瞬时速度场,在运用PIV测量特点的基础上,本文总结了以往多种湍流耗散率的测量方法,并详细介绍了大涡PIV法.对采用量纲分析方法、由耗散率定义简化而来的速度梯度方法,以及采用不同亚网格模型的大涡PIV方法的结果作了详细的比较,发现,运用不同亚网格模型的大涡PIV方法所得到的结果在槽道中心部分的结果与DNS结果非常接近,仅在壁面附近有一定的差别,这是因为在壁面附近亚网格量对耗散率的估算影响较大,不同的亚网格模型间有较大差异.采用壁面修正的Smagorinsky模型在上述方法中具有最高的精度.

普通的PIV测量得到的是一个平面上的速度

Fig.4　Distributionofturbulentdissipationrate

estimatedbydifferentmethods

分布,不能得到空间的三维速度场,而大涡PIV方法中的Smakorinsky模型或梯度模型都需要速度

场的三维的空间信息,在现有的测试手段下只能用二维的测量结果去近似处理三维速度场,这也是该方法造成误差的重要原因.因此,对于近似二维的或各向同性的流动,采用大涡PIV方法估算湍流耗散率才是可行的.

PIV由于其自身的特点,也存在一定的问题:首先,PIV查问区间限制了其空间分辨率,使得脉动速度空间导数的精确测量以及近壁面附近的测量十分困难,高分辨率的PIV算法或PTV算法是解决这一问题的关键;其次,传统的PIV只能运用于单相速度测量,这限制了其应用范围.通过图像处理的方法来扩展PIV应用范围,使其能运用于两相流动测量,在目前是一个可行的方法.

符　号　说　明

　A———Smagorinsky模型中的修正系数

k———湍流动能LCCD———CCD的像素数LIC———查问区间的像素数Lr———CCD可分辨的最小尺寸Lv———CCD所拍摄的实际范围l———湍流积分长度尺度Sij— ——应变率张量

U i, Pi———大涡模拟中所表示的宏观速度、压力ui, pi———大涡模拟中的亚网格速度、压力u,v,w———3个方向的瞬时速度

采用壁面修正的Smagorinsky模型和梯度模型

的大涡PIV方法所得结果在距壁面较远处所得结果基本重合,且与DNS结果吻合很好.只是在壁面附近,采用梯度模型的大涡PIV方法的结果出现一非常大的峰值,与DNS有较大差别,带有壁面修正的Smagorinsky模型的结果与DNS结果最为接近,但其在Y

为20时(相当于距壁面1.4

mm)就出现了明显增大的趋势,并在壁面出现约

是DNS结果两倍的峰值.除了估算方法本身的不足外,造成这一偏差的原因还有两个:首先,DNS模拟的是理想情况下的槽道流动,而实验测量是有壁面粗糙度影响的,这将影响到槽道内,特别是近壁区耗散率的值及其分布;其次,PIV测量在近壁面区域精度会有所降低,这也将影响到壁面附近耗散率的准确估算.

从上述几种不同方法所得到的耗散率分布情况看,量纲分析的方法最为简单,但该方法在定量上会存在较大的误差.依据式(3)的速度梯度方法得到的耗散率分布与采用梯度模型的大涡PIV方法的结果在最大值的数量和出现的位置都十分接近,但由于PIV空间分辨率的制约,仍限制了该方法的精度.采用梯度模型的大涡PIV方法在壁面附近会出现较大的峰值,与DNS结果有较大差别.修正后Smagorinsky模型的大涡PIV方法在上述各种方法中具有最高的精度.

5　结　论

普通的PIV测量得到的是一个平面上的速度

Fig.4　Distributionofturbulentdissipationrate

estimatedbydifferentmethods

分布,不能得到空间的三维速度场,而大涡PIV方法中的Smakorinsky模型或梯度模型都需要速度

符　号　说　明

　A———Smagorinsky模型中的修正系数

U i, Pi———大涡模拟中所表示的宏观速度、压力ui, pi———大涡模拟中的亚网格速度、压力u,v,w———3个方向的瞬时速度

采用壁面修正的Smagorinsky模型和梯度模型

为20时(相当于距壁面1.4

mm)就出现了明显增大的趋势,并在壁面出现约

　　u′———表示脉动强度的均方根脉动速度

uτ———壁面滑移速度Δ———查问区间的尺寸ε———湍能耗散率

η———Kolmogorov长度尺度λ———流动的Taylor微尺度ν———流体运动黏性系数τ——亚网格应力张量ij—

EnergyDissipation.Meas.Sci.andTech.,1996,7:1334—1348

2　AzadRS,KassabSZ.NewMethodofObtainingDissipation.

ExperimentsinFluids,1989,7:81—87

3　TsinoberA,KitE,DracosT.ExperimentalInvestigationofthe

FieldofVelocityGradientsinTurbulentFlows.J.FluidMech.,1992,242:169—192

4　AndreopoulosY,HonkanA.ExperimentalTechniquesforHighly

ResolvedMeasurementsofRotation,Strain,andDissipationRateTensorsinTurbulentFlows.Measurement

Technology,1996,7:1462—1476

Scienceand

References

1　ElsnerJW,ElsnerW.OntheMeasurementoftheTurbulence

5　ShengJ,MengH,FoxRO.ALargeEddyPIVMethodfor

TurbulenceDissipationRateEstimation.ChemicalEngineeringScience,2000,55:4423—4434

信息与交流

我国石化业将呈五强并存之势

日前,国家正式批准由中国昊华(集团)总公司、中国蓝星(集团)总公司两大公司组建中国化工集团公司,中国化工进出口总公司更名为中国中化集团公司,自此我国石油和化工行业将出现中国石油天然气集团公司、中国石油化工集团公司、中国海洋石油总公司、中国化工集团公司、中国中化集团公司五强并存的局面.

据介绍,五强企业各具特点,主营领域优势突出.其中中石油以采为主,采炼售一体运行;中石化以炼为主,同样实行采炼售一体化;中海油以海油为本,发展空间巨大;中化工以化为主,肩负新材料产业重任;中化以进出口为依托.

五强目前的发展定位不在国内,作为代表我国产业水平的国家队,他们的舞台在全球,国际化将是其重要发展策略.据了解,其实公司名称并不能表述其现在的内涵、定位和市场价值,更不能表述公司未来的发展方向.这些公司经营的内涵都已形成完整的价值链和相互支撑的产业群,成为主营业务领域推动中国经济可持续发展的骨干力量.他们将从国际和国内两个市场寻求经营资源和发展空间,加紧进行战略转移的大趋势.

(摘自:中国化工信息网)

范文二：采用大涡PIV方法研究搅拌槽内湍流动能耗散率

2008 年 6 月 The Chinese Journal of Process Engineering June 2008

采用大涡 PIV 方法研究搅拌槽内湍流动能耗散率 11112刘心洪～闵健～潘春妹～高正明～陈文民

(1. 北京化工大学化学工程学院，北京 100029;2. 中国中煤能源集团公司，北京 100010)

摘要:在槽径为 0.476 m 的六直叶涡轮桨搅拌槽内，采用粒子图像测速仪(PIV)对桨叶区的流场进行了实验研究，得到了桨叶区的平均流速和湍流动能(k)分布，采用大涡 PIV 方法对湍流动能耗散率(ε)分布进行了估算，计算了ε 与 k 的 o o相关系数. 结果表明大涡 PIV 方法能有效地估算ε 分布;桨叶区的射流向上倾斜，两尾涡分布于射流两侧，射流的倾 o 角和两尾涡中心间距随射流向壁面运动而变化，射流倾角先增大再减小，相位角θ =40时达到最大值 13.2，两尾涡中心间距先减小再增大，θ =20时达到最小值 0.038 7(用槽径 T 无因次化);湍流动能和湍流动能耗散率峰值均位于尾涡靠近射流的区域;湍流动能和湍流动能耗散率的平均相关系数为 0.363，射流核心区相关系数小于周边区域.

关键词:粒子图像测速仪;大涡模拟;湍流动能;湍流动能耗散率;尾涡;相关系数

中图分类号:TQ027 文献标识码:A 文章编号:1009?606X(2008)03?0425?07

Rushton 桨搅拌槽内的流体流动特性已经有广泛研 1前言[3,4,7,10,11]，但还没有文献报道用大涡 PIV 方法估算槽究

搅拌反应器内湍流动能耗散率(ε ) 的大小及其分布内角度解析ε 分布. 本工作用大涡 PIV 方法估算 Rushton 是搅拌反应器内湍流结构的关键参数，它直接影响微观桨搅拌槽内 6 个角度的ε. 大涡 PIV 方法最突出的优点混合效率和反应产物分布、多相体系的介观特性(如气? 液分散体系的气泡大小及其分布、液?液分散体系的液滴大小及其分布)等，因此研究ε 的大小及其分布规律对是对测量的空间分辨率要求比直接测量低得多，由于湍搅拌反应器的优化及放大具有重要意义. 流动能耗散发生在非常小的 Kolomogorov 尺度上，现有 [1]确定湍流动能耗散有多种方法. Tsinober 等用 21 测量技术一般很难达到如此高的分辨率，从理论上讲用根热线组合探针，通过直接测量脉动速度梯度方法得到 [12][2] 直接测量方法估算ε 不合适. 根据湍流的串级理论，ε. Andreopoulos 等用 9 根热线组合探针，利用各向同性假设，将ε 的定义式简化为二维形式，再用直接测量方法求得ε. 湍流动能耗散主要在 Kolomogorov 尺度上发生，测量时达到或接近该尺度才能得到较准确的ε值. 湍流动能产生于大尺度，耗散于 Kolomogorov 尺度，在用探针直接测量ε 在空间分辨率上受限制，对流场本身也有一定的干扰. 激光- 多普勒测速仪(Laser-Doppler 介于两种尺度的惯性子区内湍流动能既不产生也不耗 Velocimeter, LDV)技术解决了测量对流场的干扰问题，散，只逐级传递，且传递的能量等于在 Kolomogorov 尺 [3,4]Wu 等用 LDV 测量了 Rushton 桨搅拌槽内的流场，度上耗散掉的能量，惯性子区尺度比 Kolomogorov 尺度

大得多. 大涡 PIV 方法正是依据这种理论，不直接在计算了ε. 粒子图像测速仪(Particle Image

Kolomogorov 尺度上求ε，而是在惯性子区尺度上求得传 Velocimeter,递的湍流动能，从而大大降低了对测量分辨率的要求. [5]PIV)技术突破了单点测量的限制，Baldi 等采用 PIV 通本工作还用二维 PIV 方法得到了每个相位角度的平均过直接测量方法计算了翼型桨搅拌槽内的ε. 量纲分析法是一种比较简单的估算ε 的方法，但只能用于粗略估流速和湍流动能 k，仔细分析了尾涡和射流，通过比较 [6][7]计. Wernersson 等用能谱法、Escudé 等用湍流动能平 k 与ε 分布，发现两者分布相似，用相关系数分析了两者[8]衡法计算了搅拌槽内的ε. Sheng 等在 2000 年提出大涡的关系. PIV 方法，在 P-4 桨搅拌槽内计算ε，得到了较好的结果. 大涡 PIV 方法与直接测量方法相比对空间分辨率要求 [9]较低. Kilander 等用大涡 PIV 方法计算了翼型桨方槽内实验装置2的ε 分布.

实验采用平底圆形有机玻璃搅拌槽，内径 T=476

mm，配置单层标准六直叶涡轮桨，直径 D=1/3T，距槽底距离 C=T/3，全挡板条件，挡板宽度为 T/10. 以清水为工作介质，液位高度 H=476 mm，如图 1 所示. 为避免圆形壁面对光的折射，将搅拌槽安置在一个同样材质

的长方形器皿内，内部充满相同液位的清水介质. 示踪

粒子为直径约 8?12 μm 的空心玻璃珠，桨叶转速 N=150

收稿日期:2008?03?10，修回日期:2008?04?16 基金项目:国家自然科学基金资助项目(编号:20776008);国家重点基础研究发展规划(973)基金资助项目(编号:2007CB714302) 作者简介:刘心洪(1982?)，男，湖北省荆州市人，博士研究生，化学工程专业;高正明，通讯联系人，Tel: 010-64418267, E-mail: gaozm@mail.buct.edu.cn.

过程工程学报第 8 卷 426

4r/min，雷诺数 Re=6.2×10. 3大涡 PIV 方法

在较高雷诺数的流动中，湍流动能产生于与积分尺

度相当的大尺度脉动，耗散于 Kolmogorov 尺度的脉动，

介于这两种尺度之间存在着范围较大的惯性子区，湍流

动能在惯性子区内既不产生也不耗散，只从大尺度向小

尺度传递. 惯性子区传递的湍流动能与在 Kolmogorov

尺度上耗散的湍流动能相等，因此，计算湍流动能耗散 [6]率只需计算出在惯性子区内传递的湍流动能.

大涡模拟用位于惯性子区内的截断尺度对 Navier? Stokes (N?S)方程进行过滤，直接计算可解尺度上的脉动，对不可解尺度上的脉动用亚格子模型求解. 把瞬时图 1 实验装置图量分解为可解尺度上的量和不可解尺度上的量: Fig.1 Scheme of the experimental apparatus

本实验采用 TSI 公司的 PIV 系统. PIV 在间隔时间 '' '' (1) u =U +u , p = p + p , i i i 非常精确的 2 个时刻，瞬间“冻结”住流场，通过互相关算法运算，得到在间隔时间内流场中大量示踪粒子的其中，u, p 是瞬时速度和压力;U, p 是滤波后可分辨i i 位移，获得流场中一个平面内多点的速度. PIV 得到的的速度和压力，即可解尺度上的速度和压力，PIV 实验速度矢量是查问区内速度的空间平均值，查问区的边长 '' '' 得到的速度看作滤波后可分辨的速度;u，p是过滤掉 i i 代表 PIV 可分辨的最小长度. 本次实验的查问区为正方

的速度和压力，即不可解尺度上的速度和压力.形，边长为 0.648 mm.

对不可压缩 N?S 方程进行过滤，得到大涡控制方 PIV 测量区域为一长方形平面，选择迪卡尔坐标程: 系，坐标轴 r, z, ? 分别代表径向、轴向、切向，桨叶中心为原点，待测平面为轴向和径向构成的平面 . ?U i = 0,(2) Nd:Yag(钇钕石榴石)激光器发出的激光束通过柱面透 ?x i 镜形成片光，垂直于方槽的玻璃壁面射入圆柱形有机玻 2? U U ()?τi j ? U?U?p 璃槽内，片光与 z 轴重合. 用标准分辨率(4 008×2 672) ij i i ,+= ? + ν?(3) ?t ?x ?x?x ?x ?x j i j j j 的 CCD 相机拍摄流场中的粒子图像. 相机安放在与激光入射平面相邻的平面前，拍摄区域集中在桨的叶端部式(3)中，τ =u u ? U U是亚格子应力张量. 将式(3)乘 ij i j i j 及一部分射流区. 桨叶相对于光平面的夹角(相位角)为

，得到大涡动能控制方程，即可解尺度的动能输以U oooooo i θ，如图 2 所示. 在θ =0, 10, 20, 30, 40, 50六个相位运方程:角度各拍摄 400 对图像. ? ? ??k?k?k res res res ?p U+ U =? + ν ? τ U ? ? ij i j j ?? ?x ?x?t ?x j j j ?? ?U ?U i i ν +2 τ S , (4) ij ij ?x ?x j j Nd:Yag laser 式中，k是可解尺度的湍流动能，是可解尺度的变Sres ij Lightsheet 形率张量， θ ? ? ?U 1 ?U j i S .= ?+ ?(5) ij ?? ?x2 ?x j i ??

式(4)中的最后一项表示从可解尺度传递到不可解尺度 Cross-correlation camera 的动能，等于在 Kolmogorov 尺度上耗散掉的动能，因 4008 × 2672 P

此得到湍流动能耗散率ε 的计算式:

图 2 装置俯视图 (6) ε ?<ε>=? 2<τ s="">.SGSij ij Fig.2 The planform of apparatus

第 3 期刘心洪等:采用大涡 PIV 方法研究搅拌槽内湍流动能耗散率 427

只要选定合适的亚格子模型τ，通过式(6)就可求出运用局部各向同性假设，得ij

ε =9/ 5ε+ε +ε +ε +ε . ε.(12) ()11 22 33 12 21 2(7) Re=DNρ/μ, 搅拌槽内雷诺数为

4结果与讨论ρ 为工作介质水的密度，μ 为水的粘度，得到其中， 4Re=6.2×10. 搅拌槽内流体处于完全湍流状态，输入功4.1 平均速度场率表示为分别对 6 个相位角得到的速度场进行平均，得到平

3 5 P = N ρ N D,(8) 均速度 u 分布，如图 3 所示. P 4 为功率准数. Re=6.2×10时，功率准数取 4.6，其中，NP 图 3 中的 6 幅图分别是 6 个相位的平均速度场，速求得输入功率 P=7.08 W. 平均湍流动能耗散率为

度用叶端线速度 V无因次化，坐标轴用槽径 T 无因次 tip

化. 流体从桨叶的上方和下方被吸入桨叶，从桨叶前方 ε = P ρV ,(9) ( ) 喷出形成射流，射流形成于桨叶运动的后方. 从图可以 o 平均其中，V 是搅拌槽内流体体积， ε =0.083 6 W/kg.看出，θ =10时，射流形成，在桨叶上方和下方各有一个尾涡，尾涡中心速度最小，以尾涡中心代表尾涡位置，Kolmogorov 尺度

上尾涡的位置(r/T, z/T)为(0.156, 0.037 4)，下尾涡的位置 1/4 3 η = ν ε , (10) ()为(0.163, ?0.036 )，两尾涡之间的距离通过两尾涡的中心位置坐标求得，如图 4 所示，点 A 和 B 分别是 2 个尾η=0.058 8 mm. 查问区边长为 1.296 mm，是 Kolmogorov 涡的中心，AB 之距代表尾涡之间的距离，用 d 表示. 尺度的 22.0 倍，显然用直接测量方法不可取，因此采用 o 大涡 PIV 方法计算湍流动能耗散率. 选用 Smagorinsky θ =10时，d 为 0.073 7(用槽径 T 无因次化)，图 5 显示了 o o亚格子模型 d 随相位角的变化，d 在θ =10时最大，θ >10，d 迅速 o o减小，在θ =20时 d 达到最小值(0.038 7)，θ >20，d 开2 2 (11) τ = ?C?S S, ij s ij o 始增大，一直到θ =50，从图中曲线的斜率分析，在 oo 1/2 θ =10?20之间，d 随相位角的变化率最大.式中，S = (2S S ), C=0.16，为 Smagorinsky 常数，?s ij ij 2 个尾涡基本上对称分布在射流两侧，射流射出方为查问区长度，为 1.296 mm，通过二维 PIV 求得?U/?x, 11 向略微向上，与水平方向形成一定的夹角(射流倾角α)，?U/?x, ?U/?x, ?U/?x，由式(2)求出?U/?x，把这 512212233 近似等于两尾涡中心点中垂线与水平方向的夹角. 通过个速度梯度代入式(5)求得 S, S, S, S, S的值，再 1122331221 oo 2 个尾涡中心点的位置可算出α，图 5 显示了α同θ 的将这 5 个值代入式(11)和(6)求得ε, ε, ε, ε, ε. 由于 oo o 1122331221关系，在θ =10?40之间，α 都随θ 增大而增大，从斜率本工作使用的二维 PIV 系统难以求出?U/?x, ?U/?x, oo1331看出α 在θ =30?40之间增加较快，在θ =40时α 达到最大值(13.2)，θ >40，α开始下降.

?U/?x, ?U/?x，以致 S, S, S, S和ε, ε, ε, ε233213312332 13312332

[8]也难以求出. 为了得到ε，采用 Sheng 等的处理方法，

ooooooθ =0θ =10θ =20θ =30θ =40θ =50

图 3 角度解析平均速度场

Fig.3 Angle-resolved mean velocity fields

过程工程学报第 8 卷 428

14 0.075

) 0.070 o(12 A α 0.065 0.060 Jet flow 10 Jet incline direciton 0.055 angle, α Horizontal Distance, d 0.050 line 8

0.045

B 0.040 6 0.035 10 20 30 40 50

oPhase angle, θ () 图 4 尾涡射流示意图图 5 各相位角度的射流倾角和两尾涡间距

Fig.5 Jet incline angles and distances between vortexes Fig.4 Scheme of trailing vortexes and jet flow

at different phase angles

Turbulence kinetic

energy dissipation,

ooooooθ =0θ =10θ =20θ =30θ =40θ =50

图 6 角度解析湍流动能分布

Fig.6 Angle-resolved turbulent kinetic energy distributions

4.2湍流动能在而形成了高湍流动能分布 . 在 0.15<0.25, 湍流动能的表达式如下=""><0.05 区域，求出各相位湍流动能的平均值，="" 如图="" 7="" 所示.="" 可以看出，随着相位角度的增加，湍流动="" ''="" '2="" '2="" '2="" (13)="" k="1/2u" u="1" u+u="" +u="" ,()="" i="" i="" r="" z="" ?="" o="" 2能先增大后减小，在θ="30时达到最大值" 0.069="" 1(用="" v="" tip="" 其中，u',="" u',="" u'分别为径向、轴向和切向的脉动速度.="" 由无因次化).="" rz?="" 于本实验采用的是二维="" piv="" 系统，不能得到切向速度，="" 0.070="" 运用局部各向同性假设，得到切向速度为="" 0.065="" ε7="" '2="" '2="" '2="" u="1/2" u+u="" .(14)="" ?r="" z="" ()="">

0.060

将上式代入式(13)，得到用径向和轴向脉动速度计算湍 0.055 流动能的表达式: 6

0.050 '2 '2 (15) k =3 /4u +u .() r z Turbulent kinetic energy, k 0.045 Jet incline angle, 2 5 图 6 显示了 6 个相位角度用 V无因次化后的湍流 tip 10 20 30 40 50 2 动能(k/V)分布，从图可以看出，高的湍流动能分布于 o tip Phase angle, θ ( ) 尾涡靠近射流区域. 这部分流体与射流接近，具有较大图 7 各相位角度射流处平均湍流动能和平均湍流动能耗散的速度，又在尾涡的作用下具有较大的剪切变形速率， Fig.7 Average turbulent kinetic energy and average turbulence 高速流体在剪切作用下形成了各种尺度的脉动. 湍流动 kinetic energy dissipation in jet flow at different phase 能产生于大尺度脉动，由于该部分大尺度脉动的大量存 angles

第 3 期刘心洪等:采用大涡 PIV 方法研究搅拌槽内湍流动能耗散率 429

4.3 湍流动能耗散率<0.25><0.05 区域，求度的变化，在="" 0.15="" 通过式(12)求得湍流动能耗散率的分布，进行无因="" 出各相位湍流动能耗散率的平均值.="" 图="" 7="" 显示了射流区="" 23o="" o湍流动能耗散率平均值随相位角度的变化，随着θ="" 增大，="" 次化[ε/(dn)]，如图="" 8="" 所示.="" θ="10时，高湍流动能耗散湍流动能耗散率先增大后减小." 在θ="30，湍流动能耗散" 23="" 率最大，为="" 7.28(用="" dn无因次化).="" 湍流动能和湍流动="" 集中在桨叶上下="" 2="" 个尾涡区和射流区，3="" 个区域基本互="" 能耗散率随着射流向壁面运动都存在着先增大后减小="" o="" 不连通.="" 随着θ="" 增大，射流向槽壁运动，两尾涡之间距="" 离减小，θ="20时，3" 个区域合到一起，形成一个高湍="" o="" 流动能耗散率区，与湍流动能分布相似，峰值都集中在="" 的规律，且都在θ="30达到最大值." 尾涡靠近射流区域.="" 为了分析湍流动能耗散率随相位角="">

ooooooθ =0θ =10θ =20θ =30θ =40θ =50

图 8 角度解析湍流动能耗散分布

Fig.8 Angle-resolved turbulent kinetic energy dissipation distributions

18 (b) r/T=0.22 (a) r/T=0.2 14 Ref.[7]16 Ref.[10] 12 14 This work 10 12 ) ) 10 228 DD 338 6 /(N/(Nε ε6 4 4 Ref.[7] 2 2 Ref.[10] This work 0 0

0.00 0.02 -0.04 -0.02 0.04 0.06 -0.04 -0.02 0.00 0.02 0.04 0.06 z/T z/T

图 9 不同研究者所得轴向ε无因次化后的比

较 Fig 9 Comparison of axial profiles of normalized ε of other researchers 个相位角度的湍流动能耗散率平均，取Kolmogorov 尺度，导致实验结果偏小. 本实验所采用的将上述 6

[7][10]大涡 PIV 方法避开了空间分辨率的限制，但也存在一些 =0.2 和 0.22 处的值与 Escudé 等和 Baldi 等所得r/T 结果进行比较，如图 9 所示. 3 组曲线的趋势基本一致，缺陷. 首先，本实验采用的是二维 PIV 系统，运用了各

峰值均位于射流区域附近，因此本实验采用的大涡 PIV 向同性假设;其次，Smagorinsky 亚格子模型中给定的方法是可行的. 然而，湍流动能耗散率的测量是一项复 C是一个常数，而不是根据当地的流场情况确定，这都 s [7]杂的任务，任何一种方法都有局限性. Escudé 等用湍会给实验结果带来误差. 流动能平衡法计算了湍流动能耗散率，在计算过程中忽 [10] 4.4 湍流动能耗散率和湍流动能的比较略了压力速度相关项，引入一定的误差. Baldi 等用直比较湍流动能分布和湍流动能耗散率分布可以发接测量方法得到湍流动能耗散率，要求 PIV 的分辨率接

2 23近 Kolmogorov 尺度，但实验的 PIV 空间分辨率大于现，在高分布区[k/V>0.2, ε/(DN)>20]，两者的分布比 tip

较相似，都分布于尾涡靠近射流的区域. 这是由于这部

过程工程学报第 8 卷 430

分流体既具有较高的速度又具有较大的剪切形变速率，下结论:

高速流体被剪切为一系列湍涡. 大尺度湍涡产生了高的 (1) 六直叶涡轮桨桨叶扫过处产生射流，射流方向 o oα)，随射流向壁面运与水平方向形成一定夹角(射流倾角湍流动能，小尺度湍涡产生了高的湍流动能耗散率. 因动发生变化，相位角θ =40时α达到最大值 13.2，射流的上下各形成一个尾涡，两尾涡距离随射流运动而变此该区域既有较高的湍流动能又有较高的湍流动能耗 o 化，θ =20时两尾涡距离最小，为 0.038 7(用槽径 T 无因散率. 次化).为了了解湍流动能和湍流动能耗散率之间的关系，

(2) 湍流动能和湍流动能耗散率的峰值均位于尾涡求出各点湍流动能耗散率和湍流动能的相关系数，取值靠近射流的区域. 一般介于?1?1 之间，绝对值越大表示两者的关系越密(3) 湍流动能和湍流动能耗散率随着相位角的增加切. 各点湍流动能耗散率和湍流动能的相关系数计算式 o 先增大后减小，且两者都在θ =30达到最大值.如下: (4) 用相关系数计算了湍流动能和湍流动能耗散率

的关系，平均相关系数为 0.363. 在射流区相关系数普遍 400 小于其他区域，这是因为在射流区，湍动程度大于其他 (k?k )(ε ?ε )? n n n =1区域，该区域的湍流动能谱和耗散谱分离程度更大，两 R(k ,ε ) =(16) , 400 400 2 2种尺度之间的联系较其他区域更小. (k ?k )(ε ?ε )? ? n n n =1n =1 符号表:

求得平均相关系数为 0.363，如图 10 所示. 从图可以看

到，在射流区相关系数普遍小于其他区域. 由于湍流动

能产生于大尺度脉动，耗散于小尺度脉动，湍流动能主

C 桨叶离底距离 (mm) Smagorinsky 常数 Cs 要受大尺度脉动的影响，而湍流动能耗散主要由小尺度 d D 桨叶直径 (mm) 两尾涡之间的距离脉动控制，本实验的 Re 较大，湍流动能谱和耗散谱分 2 2 H k 湍流动能 (m /s ) 液位高度 (mm) 22离程度较大，两种尺度的脉动之间联系较小，因此相关 N k/s) 可解尺度湍流动能 (m桨叶转速 (r/min) res p N瞬时压力 (Pa) 可解P 功率准数不可解尺度系数较小. 在射流区，湍动程度大于其他区域，该区域 p'' p i压力 (Pa) 尺度压力 (Pa) 的湍流动能谱和耗散谱分离程度更大，两种尺度之间的 P 相关系数输入功率 (W) R ?1 搅拌槽内径 (mm) 瞬时联系较其他区域更小，所以射流区的相关系数普遍小于 T S可解尺度的变形率张量 (s ) ij速度 (m/s) 不可解尺度ui 平均速度 (m/s) u 其他区域. 速度 (m/s) u'' i脉动速度 (m/s) 可解' ui3 V 尺度速度 (m/s) 叶端水的体积 (m ) U ioαV线速度 (m/s) 射流倾角 () tip 23ηε湍流动能耗散率 (m/s) Kolmogorov 尺度 (mm) 2νμ动力粘度 (m/s) 粘度 (Pa?s) o3ρθ水的密度 (kg/m) 相位角 () 22τ亚格子应力 (m/s) ij

参考文献:

[1] Tsinober A, Kit E, Dracos T. Experimental Investigation of the Field

of Velocity Gradients in Turbulent Flows [J]. J. Fluid Mech., 1992,

242: 169?192. [2] Andreopoulos Y, Honkan A. Experiments of Rotation, Strain, and

Dissipation Rate Tensors in Turbulent Flows [J]. Measurement

Science and Technology, 1996, 7: 1462?1476. [3] Wu H, Patterson G K. Laser Doppler Measurements of Turbulence

Flow Parameters in a Stirred Mixer [J]. Chem. Eng. Sci., 1989, 44: o 2207?2221.图 10 θ =20时湍流动能与湍流动能耗散的相关系数分布 [4] Wu H, Patterson G K, Van Doorn M. Distribution of Turbulence Fig.10 Correlation coefficient distribution of turbulence kinetic oEnergy Dissipation Rates in a Rushton Turbine Stirred Mixer [J]. energy dissipation and turbulent kinetic energy at θ =20 Experiments in Fluids, 1989, 8: 153?160. [5] Baldi S, Yianneskis M. On the Direct Measurement of Turbulence 结论5 Energy Dissipation in Stirred Vessels with PIV [J]. Ind. Eng. Chem. 采用 PIV 方法对六直叶涡轮桨搅拌槽内的流场进 Res., 2003, 42: 7006?7016. [6] Wernersson E S W, Tragardh C. Scale-up of Rushton Turbine Agitated 行了测量，得到了速度分布和湍流动能分布. 由于 PIV Tanks [J]. Chem. Eng. Sci., 1999, 54: 4245?4256.空间分辨率的限制，用大涡 PIV 方法估算了湍流动能耗 [7] Escudé R, Liné A. Experimental Analysis of Hydrodynamics in a 散率，得到了较好的结果. 为了比较湍流动能耗散率和 Radially Agitated Tank [J]. AIChE J., 2003, 49(3): 585?603.湍流动能的关系，求得了两者的相关系数. 可以得到以 [8] Sheng J, Meng H, Fox R O. A Large Eddy PIV Method for Turbulence

第 3 期刘心洪等:采用大涡 PIV 方法研究搅拌槽内湍流动能耗散率 431

Dissipation Rate Estimation [J]. Chem. Eng. Sci., 2000, 55: in the Impeller Stream of a Stirred Vessel from Fluctuating Velocity 4423?4434.2659?2671.Gradient Measurements [J]. Chem. Eng. Sci., 2004, 59: [9] Kilander J, Rasmuson A. Energy Dissipation and Macro Instabilities[11] Lee K C, Yianneskis M. Turbulence Properties of the Impeller in a Stirred Square Tank Investigated Using an LE PIV Approach and Stream of a Rushton Turbine [J]. AIChE J., 1998, 44(1): 13?23. LDA Measurements [J]. Chem. Eng. Sci., 2005, 60: 6844?6856.[12] 张兆顺，崔桂香，许春晓. 湍流理论与模拟 [M]. 北京:清华大 [10] Baldi S, Yianneskis M. On the Quantification of Energy Dissipation学出版社, 2005. 75?78.

Investigation of Turbulence Kinetic Energy Dissipation Rate

in a Stirred Tank Using Large Eddy PIV Approach

11112LIU Xin-hong, MIN Jian, PAN Chun-mei, GAO Zheng-ming, CHEN Wen-min

(1. College of Chemical Engineering, Beijing University of Chemical Technology, Beijing 100029, China;

2. China National Coal Group Corporation, Beijing 100010, China)

Abstract: The flow field of impeller region in a Rushton impeller stirred tank with a diameter of 0.476 m was measured using particle image velocimeter (PIV). Average velocity and turbulent kinetic energy k distribution were obtained. Turbulence kinetic energy dissipation rate ε was estimated by using a large eddy PIV approach. Cross-correlation coefficient of ε and k was calculated to analyze the relation between them. The results show that the large eddy PIV approach could effectively estimate ε distribution. Impeller stream inclines slightly upwards, accompanied with two vortices on both sides of the stream. The incline angle of impeller stream and the distance between the centers of two vortices changes as the stream moves towards the tank wall. The incline angle increases before the oo ooophase angle θ =40, then decreases with the maximal value 13.2at θ =40. The distance between the centers of two vortices decreases before θ =20, then increases with the minimum value 0.038 7 (normalized by the tank diameter T) at θ =20. The peaks of turbulent kinetic energy and turbulence kinetic energy dissipation are both located in the vortex region near the jet flow. Average cross-correlation coefficient of ε and k is equal to 0.363, the cross-correlation coefficient in the center impeller stream is larger than that in the other parts of the region.

Key words: particle image velocimeter; large eddy simulation; turbulent kinetic energy; turbulence kinetic energy dissipation rate; trailing vortex; cross-correlation coefficient

范文三：[湍动能和湍流耗散率]湍动：湍动

[湍动能和湍流耗散率]湍动:湍动篇一 : 湍动:湍动-相关条目

当流体的Re超过一定数值时，流体就处于湍流状态，这时管道内的每一个流体质点都有是在做不规则的、杂乱无章的运动，就称为湍动。

湍_湍动 -相关条目

等离子体湍动加速等离子体湍动加速

plasma turbulence acceleration

李晓卿

等离子体的1个最重要特性是不稳定性。微小的扰动就能在等离子体中激起各种等离子体波。这种等离子体的激发态通常称为等离子体湍动。湍动元和荷电粒子碰撞会引起它们之间的能量交换,从而导致粒子加速,这种现象称为等离子体湍动加速。这种加速效应带有统计性质,和经典的费密加速类似。业已证明,等离子体激元和荷电粒子间的碰撞总是导致粒子的加速率正比于L ,L 是两激元之间的平均距离,也就是两湍动元之间的平均尺度。这种关系是普遍的,并不取决于具体的加速机制。因而湍动元尺度越小,加速效率就越高。在等离子体中,存在各种高频等离子体波,它们的波长是短的,所以,加速效率就比费密加速效率大得多。计算表明,如果太阳缓变射电是由等离子体中的电子振汤波转化来的,那么,这种电子波就能在

一天之内把足够多的粒子加速到具有相当于1个耀斑爆发的能量。可见,这种湍动加速效率是非常高的。等离子体湍动加速通常包括2种情况:如果等离子体波的相速度大于粒子的热运动速度,那么,这种等离子体波只能加速少数快粒子,这叫作等离子体纯粹加速,如果波的相速度小于被加速粒子的热运动速度,那么,大多数粒子都能被这种等离子体波加速,这叫作等离子体湍动加热。

对于活动星系核)类星体)脉冲星)蟹状星云等,不管它们的辐射机制如何,为了得到和观测资料一致的结果,总得假定它们的高能粒子数随能量的分布是采取幂指数形式的。正是考虑到等离子体湍动加速效应,才有可能自洽地获得粒子的这种幂律谱。

篇二 : 第4章湍流动能稳定度和尺度之2

第四章湍流动能、稳定度和尺度

LOGO

第四章湍流动能、稳定度和尺度

1 湍流动能TKE

稳定度

尺度

Company Logo

4.1 湍流动能TKE收支方程

1 2 2 2 e ? u? ? v? ? w? 2

? ?

?e g ?u ? w‘e 1 ? w‘p‘ ? w‘,v‘ ? u‘w‘ ? ? ?ε ?t ,v ?z ?z ρ ?z

? ? ? ? ?

? ?

Company Logo

4.2 大气稳定度

4.2.1 静力稳定度 4.2.2 动力稳定度

4.2.3 湍流稳定度的判定

Company Logo

4.2 大气稳定度

4.2.1 静力稳定度

静止大气中产生浮力对流的潜在能力与气流运动无关，仅与空气密度有关

非局部定义:分析整个边界层的稳定度

Company Logo

4.2.2 动力稳定度

与气流运动有关，即使空气处于静力稳定，风切变也能产生动力湍流。

当风切变超过一定临界值，做层流运动的流体会变得动力不稳定。

开尔文,赫姆霍兹不稳定波

参考Stull书P179,181

Company Logo

开尔文,赫姆霍兹不稳定波

密度小密度大

Company Logo

澳大利亚杰维斯湾碎浪花云上方暖气流，速度快下方冷气流，速

度慢

美国落基山碎浪花云具有纤维组织的卷云像羽毛、发丝、马尾等

Company Logo

土星

Mount Duval 澳大利亚

San Francisco 美国

Company Logo

?如何理解:

静力不稳定气流总是动力不稳定的;而静力稳定时，大气由于切变运动可产生动力不稳定。

流体起着使不稳定破坏的作用。化学上的勒夏特利埃定律―如果把某一应力作用于处于平衡状态的系统，该系统将向尽量减少应力影响的方向上变化。‖

静力不稳定:对流使浮力流体向上移动，从而使系统稳定下来动力不稳定:湍流会使风切变减弱，从而也能使系统稳定下来

湍流起着消灭自身的作用。

Company Logo

4.2.3 湍流稳定度的判定

而浮力作功项则可正可负，为什么,

?u 大气边界层内雷诺应力作功项 ? u?w? 始终为正值， ?z

上升运动:

u? ? u?z ? ? u?z ? dz ? ? 0

? u?w? ? 0

下沉运动:

u? ? u?z ? ? u?z ? dz ? ? 0

? u?w? ? 0

动量输送始终为正值～

Company Logo

白天，剪切产生的湍流尺度偏小，浮力产生的湍流尺度偏大，对 / 错,

错～剪切贡献于水平湍流通量，而水平方向无限制;

浮力贡献于垂直方向湍流通量，而垂直方向下有下垫面、上有ABL顶的限制;

因此:湍流由剪切和浮力两项共同作用，始终由平均正向传输给湍流能量，非各向同性

引出稳定度判据

Company Logo

?u ? u j ? e 1 ? ui ? p‘ ?e ?e g ?uj ? ui ?? ?? i 3 ? ui ? u j ? i ? ? ?? ?t ?x

j ? ?x j ?x j ? ?xi

? ?

取第湍流能量方程

?、?两项的比值可以得到一个重要的判据:

1) 通量理查孙数 Rf

g w‘,v‘ ,v 浮力项 Rf ? ~ ?ui 雷诺应力项 ui‘u j‘ ?x j

水平均匀无下沉

g w‘,v‘ ,v ?u ?v u‘w‘ ? v‘w‘ ?z ?z

Rf 0:静力稳定。

Company Logo

理论上: Rf 1:动力稳定，气流变成片流。

实际上:

当湍流在较小的Rf 时即不能维持，通常取临界值 0.25

Company Logo

2) 梯度理查逊数 Ri

?? K 理论: w‘? ? ? ? K h ?z

g w‘,v‘ ,v Rf ? ?u ?v u‘w‘ ? v‘w‘ ?z ?z

?u w‘ u ? ? ? K m ?z

假设:Km = Kh

g ?? v ? v ?z ? ?u ? ? ?v ? ? ? ? ? ? ?z ? ? ? ?z ? ? ? ? ?

2 2

0:静力稳定

Company Logo

Ri ?

用梯度理查逊数Ri判定静力、动力稳定度

令湍流开始时的 Ri 临界值为 Rc ，湍流终止时的 Ri 临界值为 RT

Ri RT ，动力稳定，湍流变成层流，一般RT=1.0

Ri 静力稳定度动力稳定度气流状态 -2 -1 0 0.25 稳定依赖于过去稳定层流 1 2

不稳定不稳定湍流

依赖于过去

Company Logo

Ri ?

g ?? v ? v ?z ? ?u ? ? ?v ? ? ? ?? ? ? ?z ? ? ?z ? ? ? ? ?

2 2

3) 总体理查逊数 RB:

g ?? v ? v ?z ? ?u ? ? ? v ? ? ? ?? ? ? ?z ? ? ?z ? ? ? ? ?

2 2

RB ?

具有平滑作用，临界值大于0.25, ?Z越薄，临界值越接近于0.25

Company Logo

100%

晴空湍流

0 0 0.25 10.25

总体理查逊数注意:理查逊数只能表示有或无湍流、并不能说明湍流强度。

Company Logo

? 个例分析——夜间Ri 的演变

阴影区为Ri * ?

m *

4.3.4 其它尺度参数

Company Logo

?ui ? u j ? e 1 ? ui ? p‘ ?e ?e g ?uj ? ui ?? ?? i 3 ? ui ? u j ? ? ? ?? ?t ?x

j ? ?x j ?x j ? ?xi

? ?

4.3.1 Monin,Obukhov长度

M-O长度是近地层的尺度参数

TKE方程两边同乘以

...... ? ? k zg w‘? v ? s

k卡门常数 0.35-0.42

3 v *

k z ui ? u j ? s ?ui k z? s ? ? ... ? 3 3 ?x j u* u*

Company Logo

稳定度参数,

z ? k zg w‘? v ? s ? ? ? 3 L ? v u*

?? 0 静力稳定 ? ?? 0 静力中性 ?? 0 静力不稳定 ?

M-O 长度

? ? v u* L? k g w‘? v ? s

M-O长度:与地面以上某一高度成正比，在这个高度上，浮力因子首先超过湍流的机械产生。

Company Logo

?u L?? ?g w‘? v ?

3 v *

?? 0 不稳定 z ? ? ? ? ?? 0 中性 L ? ?? 0 稳定

z L

>>1 浮力占优势 = 1 浮力与雷诺应力作功相等的物理意义及变化

Company Logo

5、大气稳定层结的判据: 通量Richardson数、梯度Richardson数、整体Richardson数由Ri判定大气稳定度 Monin,Obukhov长度

6、闭合理论:

闭合理论是什么,为什么出现闭合问题,

Company Logo

?下一章内容:定常条件下的大气边界层 ?前提:复习流体力学中相似原理与量纲分析一章的内容

Company Logo

LOGO

篇三 : 介绍湍流 54

第三章，湍流模型

第一节，前言

湍流流动模型很多，但大致可以归纳为以下三类:

第一类是湍流输运系数模型，是Boussinesq于1877年针对二维流动提出的，将速度脉动的二阶关联量表示成平均速度梯度与湍流粘性系数的乘积。即:

u2=μt?ρu1

?u1

3,1 ?x2

推广到三维问题，若用笛卡儿张量表示，即有:

??ρkδij 3,2 ?3?

模型的任务就是给出计算湍流粘性系数μt的方法。根据建立模型所需要的微分方程的数

目，可以分为零方程模型，单方程模型和双方程模型。

第二类是抛弃了湍流输运系数的概念，直接建立湍流应力和其它二阶关联量的输运方程。第三类是大涡模拟。前两类是以湍流的统计结构为基础，对所有涡旋进行统计平均。大涡模拟把湍流分成大尺度湍流和小尺度湍流，通过求解三维经过修正的Navier-Stokes方程，得到大涡旋的运动特性，而对小涡旋运动还采用上述的模型。

实际求解中，选用什么模型要根据具体问题的特点来决定。选择的一般原则是精度要高，应用简单，节省计算时间，同时也具有通用性。

FLUENT提供的湍流模型包括:单方程模型、双方程模型及雷诺应力模型和大涡模拟。

RANS-based models

??ui?uj

′′+?ρuiuj=μt?

??x

?j?xi

包含更多

物理机理每次迭代计算量增加

提的模型选

湍流模型种类示意图

第二节，平均量输运方程

雷诺平均就是把Navier-Stokes方程中的瞬时变量分解成平均量和脉动量两部分。对于速度，有:

ui=ui+ui′ 3,3

其中，ui和ui′分别是平均速度和脉动速度

类似地，对于压力等其它标量，我们也有:

θ=θ+θ′ 3,4 其中，θ表示标量，如压力、能量、组分浓度等。

把上面的表达式代入瞬时的连续与动量方程，并取平均，我们可以把连续与动量方程写成如下的笛卡儿坐标系下的张量形式:

?ρ?

+=0 3,5 ?t?xi

Dui?p?ρ=?+Dt?xi?xj

???ui?uj2?ul???

??+?ρuiuj 3,6 +?δij?μ???????xj?xi3?xl????xj

上面两个方程称为雷诺平均的Navier-Stokes方程。他们和瞬时Navier-Stokes方程有相同的形式，只是速度或其它求解变量变成了时间平均量。额外多出来的项?ρuiuj是雷诺应力，表示湍流的影响。如果要求解该方程，必须模拟该项以封闭方程。

如果密度是变化的流动过程如燃烧问题，我们可以用法夫雷平均。这样才可以求解有密度变化的流动问题。法夫雷平均就是出了压力和密度本身以外，所有变量都用密度加权平均。变量的密度加权平均定义为:

Φ=ρΦ/ρ 3,7 符号,表示密度加权平均;对应于密度加权平均值的脉动值用Φ′′表示，即有:~

Φ=Φ+Φ′′。很显然，这种脉动值的简单平均值不为零，但它的密度加权平均值等于零，即:

Φ′′?0， ρΦ′′=0

Boussinesq近似与雷诺应力输运模型

为了封闭方程，必须对额外项雷诺应力?ρuiuj进行模拟。一个通常的方法是应用Boussinesq假设，认为雷诺应力与平均速度梯度成正比，即:

??ui?uj

?ρuiuj=μt?+

??x

?j?xi

??δij 3,8 ?3?xi?

Boussinesq假设被用于Spalart-Allmaras单方程模型和k?ε双方程模型。Boussinesq近似

的好处是与求解湍流粘性系数有关的计算时间比较少，例如在Spalart-Allmaras单方程模型中，只多求解一个表示湍流粘性的输运方程;在k?ε双方程模型中，只需多求解湍动能k和耗散

率ε两个方程，湍流粘性系数用湍动能k和耗散率ε的函数。Boussinesq假设的缺点是认为湍流粘性系数μt是各向同性标量，对一些复杂流动该条件并不是严格成立，所以具有其应用限制性。

另外的方法是求解雷诺应力各分量的输运方程。这也需要额外再求解一个标量方程，通常是耗散率ε方程。这就意味着对于二维湍流流动问题，需要多求解4个输运方程，而三维湍流问题需要多求解7个方程，需要比较多的计算时间，对计算机内存也有更高要求。

在许多问题中，Boussinesq近似方法可以得到比较好的结果，并不一定需要花费很多时间来求解雷诺应力各分量的输运方程。但是，如果湍流场各向异性很明显，如强旋流动以及应力驱动的二次流等流动中，求解雷诺应力分量输运方程无疑可以得到更好的结果。

第三节，湍流模型

3.3.1 单方程模型

~，表征出了近壁区域以外的湍流运动粘Spalart-Allmaras模型的求解变量是ν

~的输运方程为: 性系数。ν

~~?~????νDν1????ν??~???Yν 3,9 ρ=Gν+???+Cb2ρ???Dtζν~??xj?????xj????xj??

其中，Gν是湍流粘性产生项;Yν是由于壁面阻挡与粘性阻尼引起的湍流粘性的减少;ζν~

和Cb2是常数;ν是分子运动粘性系数。

湍流粘性系数用如下公式计算:

~f μt=ρνν1

其中，fν1是粘性阻尼函数，定义为:fν1湍流粘性产生项，Gν用如下公式模拟:

~χ3ν

=3，并且χ?。

νχ+Cν31

Gν=Cb1ρSν 3,10

~νχ~

其中，S?S+22fν2，而fν2=1?。其中，Cb1和k是常数，d是计算点到

1+χfν1kd

2?ij?ij。?ij定义为:

? 3,11 ??

?uj?ui1?

?ij=??

?2??xi?xj

壁面的距离;S?

由于平均应变率对湍流产生也起到很大作用，FLUENT处理过程中，定义S为:

S??ij+Cprodmin 3,12

其中，Cprod=2.0，?ij?

?ij?ij，Sij?2SijSij，平均应变率Sij定义为:

Sij=

?u1??j+?ui2???xi?xj?

? 3,13 ??

在涡量超过应变率的计算区域计算出来的涡旋粘性系数变小。这适合涡流靠近涡旋中心的区域，那里只有―单纯‖的旋转，湍流受到抑止。包含应变张量的影响更能体现旋转对湍流的影响。忽略了平均应变，估计的涡旋粘性系数产生项偏高。

湍流粘性系数减少项Yν为:

~??ν

Yν=Cw1ρfw?? 3,14

?d?

?1+Cw?3

3,15 其中，fw=g?6

?g+Cw3?

g=r+Cw2 3,16

~ν

r?22 3,17

Skd

~ν~

其中，Cw1，Cw2，Cw3是常数，S?S+22fν2。在上式中，包括了平均应变率对S

的影响，因而也影响用S计算出来的r。

上面的模型常数在FLUENT中默认值为:Cb1=0.1335，Cb2=0.622，ζν~=2/3，

1/6

Cν1=7.1，Cw1=Cb1/k2+/ζν~，Cw2=0.3，Cw3=2.0，k=0.41。

壁面条件

~设置为零。当计算网格足够细，可以计算层流底层时，壁面

切应在壁面，湍流运动粘性ν

力用层流应力,应变关系求解，即:

ρuyu

=η 3,18 uημ

如果网格粗错不能用来求解层流底层，则假设与壁面近邻的网格质心落在边界层的对数区，则根据壁面法则:

?ρuηy?u1

=lnE?? 3,19 ?μ?uηk??

其中，k=0.419，E=9.793。

对流传热传质模型

在FLUENT中，用雷诺相似湍流输运的概念来模拟热输运过程。给出的能量方程为:

cpμt?????

+[ui]=???k+Prt?t?xi?xi???

???T

???x+ujeff?+Sh 3,20

??i?

式中，E是总能量，eff是偏应力张量，定义为:

eff=μeffeff表示粘性加热，耦合求解。如果默认为分开求解，FLUENT不求解处eff。但是可以通过变化―粘性模型‖面板上的湍流普朗特数，其默认值为0.85。

湍流质量输运与热输运类似，默认的Schmidt数是0.7，该值同

样也可以在―粘性模型‖面板上调节。

标量的壁面处理与动量壁面处理类似，分别选用合适的壁面法则。

综上所述，Spalart-Allmaras模型是相对简单的单方程模型，只需求解湍流粘性的输运方程，并不需要求解当地剪切层厚度的长度尺度。该模型对于求解有壁面影响流动及有逆压力梯度的边界层问题有很好模拟效果，在透平机械湍流模拟方面也有较好结果。

Spalart-Allmaras模型的初始形式属于对低雷诺数湍流模型，这必须很好解决边界层的粘性影响区求解问题。在FLUENT中，当网格不是很细时，采用壁面函数来解决这一问题。当网格比较粗糙时，网格不满足精确的湍流计算要求，用壁面函数也许是最好的解决方案。另外，该模型中的输运变量在近壁处的梯度要比k?ε中的小，这使得该模型对网格粗糙带来数值误差不太敏感。

但是，Spalart-Allmaras模型不能预测均匀各向同性湍流的耗散。并且，单方程模型没有考虑长度尺度的变化，这对一些流动尺度变换比较大的流动问题不太适合。比如，平板射流问题，从有壁面影响流动突然变化到自由剪切流，流场尺度变化明显。

3.3.2 标准k?ε模型

标准k?ε模型需要求解湍动能及其耗散率方程。湍动能输运方程是通过精确的方程推导得到，但耗散率方程是通过物理推理，数学上模拟相似原形方程得到的。该模型假设流动为完全湍流，分子粘性的影响可以忽略。因此，标准k?ε模型只适合完全湍流的流动过程模

拟。

标准k?ε模型的湍动能k和耗散率ε方程为如下形式:

??μt??k???μ+???+Gk+Gb?ρε?YM 3,22 ?ζk??xi???

μt??ε?Dεεε2???

? 3,23 ρ=μ+???+C1ε?C2ερ??Dt?xi??ζk??xi?kk

在上述方程中，Gk表示由于平均速度梯度引起的湍动能产生，Gb是用于浮力影响引起Dk?

=Dt?xi

的湍动能产生;YM可压速湍流脉动膨胀对总的耗散率的影响。湍流粘性系数μt=ρCμ。

在FLUENT中，作为默认值常数，C1ε,1.44，C2ε=1.92，Cμ=0.09，湍动能k与耗散

率ε的湍流普朗特数分别为ζk,1.0，ζε,1.3。可以通过调节―粘性模型‖面板来调节这些常数值。

介绍湍流 54_湍流动能

3.3.3 重整化群κ-ε模型

重整化群κ-ε模型是对瞬时的Navier-Stokes方程用重整化群的数学方法推导出来的模型。模型中的常数与标准κ-ε模型不同，而且方程中也出现了新的函数或者项。其湍动能与耗散率方程与标准κ-ε

模型有相似的形式:

??k?αμ?keff?+Gk+Gb?ρε?YM 3,24 ?xi??

Dε???ε?εε2

ρ=?R 3,25 ??+C1ε?C2ερDt?xi??xi?kk

Gk表示由于平均速度梯度引起的湍动能产生，Gb是用于浮力影响引起的湍动能产生;YM可压速湍流脉动膨胀对总的耗散率的影响，这些参数与标准κ-ε模型中相同。αk和αε分别是湍动能k和耗散率ε的有效湍流普朗特数的倒数。 ρDk?=Dt?xi

湍流粘性系数计算公式为:

~?ρ2k?ν~ 3,26 ??d=1.72dν??3?1?Cν??~=μ/μ，C?100 其中，νeffν

对上面方程积分，可以精确得到有效雷诺数对湍流输运的影响，这有助于处理低雷诺数和近壁流动问题的模拟。

对于高雷诺数，上面方程可以给出:μt=ρCμ，Cμ=0.0845。这个结果非常有意ε

思，和标准κ-ε模型的半经验推导给出的常数Cμ=0.09非常近似。

在FLUENT中，如果是默认设置，用重整化群κ-ε模型时候是针对的高雷诺数流动问题。如果对低雷诺数问题进行数值模拟，必须进行相应的设置。

重整化群κ-ε模型有旋修正

通常，平均运动有旋时候对湍流有重要影响。FLUENT中重整化群κ-ε模型通过修正湍流粘性系数来考虑了这类影响。

湍流粘性的修正形式为:

kμt=μt0f 3,27 ε

其中，μt0是不考虑有旋计算出来的湍流粘性系数;Ω是FLUENT计算出来的特征旋流数;αs是旋流常数，不同值表示有旋流动的强度不同。流动可以是强旋或者中等旋度的。FLUENT默认设置αs,0.05,针对中等旋度的流动问题，对于强旋流动，可以选择较大的值。

湍动能及其耗散率的有效湍流普朗特数倒数的计算公式为:

α?1.3929

α0?1.39290.6321α+2.3929α0+2.39290.3679=μmol 3,28 μeff

式中，α0,1，在高雷诺数流动问题中，μmol/μeff??1，αk=αε=1.393。

湍流耗散率方程右边的R为:

k1+β,3

其中，,?Sk/ε，,0=4.38，β=0.012。

Dε?ρ=Dt?xi

则:C2ε=C2ε+*R=Cμρ,3ε2 3,29 为了更清楚体现R对耗散率的影响，我们把耗散率输运方程重写为: ??ε?εε2ε2*?C2ερ 3,30 ??+C1ε?C2ερ?xi?kkk?Cμρ,3

1+β,3

** 3,31 在, 和标准κ-ε模型中给出的C2ε,1.92接近。因此，对于弱旋和中等旋度的流动问题，重整化群κ-ε模型给出的结果比标准κ-ε模型的结果要大。

重整化群模型中，C1ε=1.42，C2ε=1.68。

3.3.4 可实现κ-ε模型

可实现κ-ε模型的湍动能及其耗散率输运方程为:

??μt??k????+Gk+Gb?ρε?YM 3,32 ??μ+?ζk??xj?????

μt??ε?Dε???ε2ε?ρμ+CS?C+CC3εGb 3,33 =+ρερ???121ε??Dt?xj?ζt??xj?kk+???Dk?=ρDt?xj

其中，C1=max?0.43,?

?,?，,=Sk/ε ,+5??

在上述方程中，Gk表示由于平均速度梯度引起的湍动能产生，Gb是用于浮力影响引起的湍动能产生;YM可压速湍流脉动膨胀对总的耗散率的影响。C2和C1ε是常数;ζk，ζε分别是湍动能及其耗散率的湍流普朗特数。在FLUENT中，作为默认值常数，C1ε,1.44，C2=1.9，ζk,1.0，ζε,1.2。

可实现κ-ε模型的湍动能的输运方程与标准κ-ε模型和重整化群κ-ε模型

有相同的形式，只是模型参数不同。但耗散率方程有较大不同。首先耗散率产生项不包含湍动能产生项Gk，现在的形式更能体现能量在谱空间的传输。另外的特色在于耗散率减少项中，不具有奇异性。

并不象标准κ-ε模型模型那样把K放在分母上。

该模型适合的流动类型比较广泛，包括有旋均匀剪切流，自由流，腔道流动和边界层流动。对以上流动过程模拟结果都比标准κ-ε模型的结果好，特别是可实现κ-ε模型对圆口射流和平板射流模拟中，能给出较好的射流扩张角。

湍流粘性系数公式为μt=ρCμ，这和标准κ-ε模型相同。不同的是，在可实现κ-ε

ε模型中，Cμ不再是个常数，而是通过如下公式计算: U*KA0+Asε~其中，U*=SijSij+?ij?ij，?ij=?ij?2εijkωk，?ij=?ij?εijkωk，?ij是 is the mean

rate-of –rotation tensor viewed in a rotating reference frame with

the angular velocity ωk。模型常数A0=4.04，As=Cμ=1 3,34 cosθ，而:

SijSjkSkj11?uj?ui~θ=arccos，式中W,，S=SijSij，Sij=]=??+k+ueffjijeff?+Sh 3,39 ??t?xi?xi??xi?

式中，E是总的能量，keff是有效导热系数;eff是偏应力张量，定义为:

eff??uj?ui+=μeff???x?i?xj?2??μeff?uiδij 3,40 ?3?xi?

eff表示的是粘性加热，耦合求解时总是计算。如果不是耦合求解时候，作为默认设置，并不求解该量。如果有需要，需在―粘性模型‖面板中设置。

对于重整化群κ-ε模型，有效导热系数为:

keff=αcpμeff 3,41 α用计算，式中，α0=1/Pr=k/μCp。事实上，α随着μmol/μeff的变化而变化，这是重整化群κ-ε模型的一个优点，因为实验中证明，湍流普朗特数随分子普朗特数及湍流而变化。

湍流质量输运处理过程与能量输运过程类似。对于标准κ-ε模型和可实现的κ-ε模型，默认的Schmidt数是0.7，重整化群模型中，是通过方程3,28来计算的，其中，α0=1/Sc，Sc是分子Schimidt数。

3.3.5 雷诺应力模型

雷诺应力模型是求解雷诺应力张量的各个分量的输运方程。具体形式为: ??+= ?t?xk

对流项Cij

???ρuiujuk+p+?xk?xk???uiuj? ?μ?xk??

TL湍流扩撒项Dij 分子扩散Dij []

?Uj??Ui??uuuu+?ρβgiuj,+gjui, ?ρ?jk?ik?x??xk?k?

应力产生项Pij 浮力产生项目Gij

??u?uj?uj???2μ?ui+p?i+ ??x??xk?xk?j?xi?

压力应变项Φij 耗散项εij

?2ρ?kujumεikm+uiumεjkm 3,42

系统旋转产生项Fij

上面方程中，Cij，Dij，Pij，Fij不需要模拟，而Dij，Gij，Φij，εij需要模拟以封闭方程。下面简单对几个需要模拟项的模拟。

TDij可以用Delay and Harlow [L38]的梯度扩散模型来模拟，

即: LT)

?D=Cs?xkTij?kukul?uiuj??ρ? 3,43 ε?xl????

但这个模型会导致数值不稳定，因此FLUENT程序中采用标量湍流扩散模型:

?? 3,44 ??

式中，湍流粘性系数用μt=ρCμ来计算，根据Lien and Leschziner

[L98]，ζk=0.82，ε?D=?xkTij?μt?uiuj??ζk?xk?

这和标准κ-ε模型中选取1.0有所不同。

根据Gibson and Launder [L58], Fu [L55], Launder [L88,L89], 压力应变项Φij可以分解为三项，即:

wΦij=Φij,1+Φij,2+Φij 3,45

介绍湍流 54_湍流动能

Φij,1，Φij,2和Φwij分别是慢速项，快速项和壁面反射项。

2??uu?δkijij?，常数C1=1.8。 ?3??

12??Φij,2=?C2??δij?，C2=0.60，P=Pkk，23??

11G=Gkk，C=Ckk。 22Φij,1=?C1ρ

壁面反射项用于重新分布近壁的雷诺正应力分布，主要是减少垂直于壁面的雷诺正应力，增加平行于壁面的雷诺正应力。该项模拟为: εk

Φwijε?33k3/2?=C1′?ukumnknmδij?uiuknjnk?ujuknink?

k?22?Clεd

33k3/2??′?Φkm,2nknmδij?Φik,2njnk?Φjk,2nink?+C2 3,46 22Cεd??l

′,0.5，C2′=0.3，nk是xk在垂直于壁面方向上的单位分量，d是到壁面的距式中，C1

离;Cl=Cμ3/4/k，Cμ=0.09，k=0.41。

默认设置时候，FLUENT不计算Φwij。如果需要计算时候，在―粘性模型‖面板中设置。

线性压力应变模型

对于小雷诺数流动，特别是用双层模型求解近壁流动问题时，FLUENT中通过改进模型常

′和C2′来改进压力应变项Launder [L91]。这一过程只有在选择双层流模型时候，数C1，C2，C1

在―粘性模型‖面板上调节。

?2?C1=1+2.58AA2??1?e? ??

C2=0.75A

2C1′=?C1+1.67 3

1??2?C2??36′=max?C2,0? ?C2?????

2其中，Ret=ρk/，参数A和张量不变量A2, A3定义为:

?9?A??1?? ?8?

A2?αikαki

A3?αikαkjαji

式中，αij是雷诺应力张量各向异性部分，定义为:

2???ρuiuj+ρkδij3αij=??ρk?????? 3,47 ???

二阶压力应变模型

二阶压力应变模型由Spezible {L157}等人提出。

1*Φij=?bij+C2ρε+ρkSij 3

2+C4ρk+C5ρk 3,48 3

式中，bij是雷诺应力各向异性张量，定义为:

??? 3,49 ???

?uj?ui??uj?ui?1?1?????;模型常数平均应变率Sij定义为: Sij=;?ij=+????2??xi?xj?2??xi?xj??

*C1=3.4，C1*=1.8，C2=4.2，C3=0.8，C3=1.3，C4=1.25，C5=0.4。

二阶压力应变模型不需要考虑壁面反射影响去模拟对数区湍流边界层过程。

浮力对湍流的影响

浮力引起的产生项模拟为: 2???ρuiuj+ρkδij3 bij=??2ρk???

μt??T?T??? 3,50 +gjGij=βgi?xi?Prt???xj?

其中，Prt是能量的湍流普朗特数，默认设置值为0.85。

对于理想气体，把热膨胀系数的定义代入上式，得:

耗散项εij的模拟

μt??ρ?ρ??? 3,51 Gij=?gi+gj?ρPrt??xj?xi??35

耗散张量εij模拟为:

2εij=δij 3,52 3

2式中，YM=ρε2Mt，Mt是马赫数;标量耗散率ε用标准k-ε模型中的采用的耗散率输运方

程求解。

雷诺应力模型的边界条件

在流场进口，雷诺应力模型需要各个雷诺应力分量和湍动能耗散率的值。这些值可以直接输入，也可以湍流强度和特征长度来计算。

在壁面，雷诺应力模型通过壁面函数，给出各个雷诺应力分量和耗散率的值。

雷诺应力模型的能量与质量输运方程

在雷诺应力模型中，对流传热传质模型都是通过雷诺相似湍流动量输运方程得到的。能量方程形式为:

cpμt?T???[ui]=??++ujeff??t?xi?xi?Prt?xi

式中，E是总的能量;eff是偏应力张量，定义为:

eff??uj?ui+=μeff???x?i?xj???+Sh 3,53 ??2??μeff?uiδij 3,54 ?3?xi?

eff表示的是粘性加热，耦合求解时总是计算。如果不是耦合求解时候，作为默认设置，并不求解该量，并且Prt=0.85。如果有需要，需在―粘性模型‖面板中设置。

3.3.6 大涡模拟

湍流中包含了不同时间与长度尺度的涡旋。最大长度尺度通常为平均流动的特征长度尺度。最小尺度为Komogrov尺度。

LES的基本假设是1，动量、能量、质量及其它标量主要由大涡输运;2，流动的几何和边界条件决定了大涡的特性，而流动特性主要在大涡中体现;3，小尺度涡旋受几何和边界条件影响较小，并且各向同性;大涡模拟过程中，直接求解大涡，小尺度涡旋模拟，从而使得网格要求比DNS低。

3.3.6.1大涡模拟的控制方程

LES的控制方程是对Navier-Stokes方程在波数空间或者物理空间进行过滤得到的。过滤的过程是去掉比过滤宽度或者给定物理宽度小的涡旋，从而得到大涡旋的控制方程。

过滤变量定义为:

θ=?DθGdx′ 3,55 其中，D表示流体区域;G是决定涡旋大小的过滤函数。

在FLUENT中，有限控制体离散本身暗中包括了过滤运算，

θ=1?Vθdx′， x′?V 3,56 V

其中V是计算控制体体积，过滤函数为:

x′?V?1/VG=? 3,57 0′x?V?

目前，大涡模拟对不可压流动问题得到较多应用，但在可压缩问题中的应用还很少，因此这里涉及的理论都是针对不可压流动的大涡模拟方法。在FLUENT中，大涡模拟只能针对不可压流体的流动。

过滤不可压的Navier-Stokes方程后，可以得到LES控制方程:

?ρ?ρui=0 3,58 +?t?xi

?u????p?ηij+=?? 3,59 ?t?xj?xj?xj?xi?xj

其中，ηij为亚网格应力，定义为:

ηij=ρuiuj?ρui?uj 3,60 很明显，上述方程与雷诺平均方程很相似，只不过大涡模拟中的变量是过滤过的量，而非时间平均量，并且湍流应力也不同。

3.3.6.2 亚网格模型

由于LES中亚网格应力项是未知的，并且需要模拟以封闭方程。目前，采用比较多的亚网格模型为涡旋粘性模型，形式为:

1ηij?ηkkδij=?2μtSij 3,61 3

式中，μt是亚网格湍流粘性系数;Sij是求解尺度下的应变率张量，定义为:

?? 3,62 ??

求解亚网格湍流粘性系数μt时，FLUENT提供了两种方法。第一，Smagorinsky-Lilly模型;第二，基于重整化群的亚网格模型。

最基本的亚网格模型是Smagorinsky [L145]最早提出的，Lilly [L99]把它进行了改善，这就 37 ?u?uj1?Sij=?i+2???xj?xi

是今天的Smagorinsky-Lilly模型。该模型的涡粘性计算方程为:

μt=ρL2

sS 3,63

式中，Ls是亚网格的混合长度;S?2SijSij。Cs是Smagorinsky常数，则亚网格混合长度Ls可以用下式计算。

Ls=min 3,64 其中，k=0.42，d是到最近壁面的距离，V是计

算控制体体积。

Lilly通过对均匀各向同性湍流惯性子区湍流分析，得到了Cs,0.23。但是研究中发现，对于有平均剪切或者过渡流动中，该系数过高估计了大尺度涡旋的阻尼作用。因此，对于比较多的流动问题，Cs,0.1有比较好的模拟结果，该值是FLUENT的默认设置值。

我们再来看看基于重整化群思想的亚网格模型。人们用重整化群理论推导出了亚网格涡旋粘性系数[L182]，该方法得到的是亚网格有效粘性系数，μeff=μ+μt，而

μeff??μs2μeff???=μ?1+H?C??3?μ???????

1/321/3 3,65 式中，μs=2SijSij，H是Heaviside函数，

x>0?x 3,66 H=? x?0?0

V是计算控制体体积;重整化群常数Crng=0.157,而常数C=100。

对于高雷诺数流动，μeff?μt，基于重整化群理论的亚网格模型就与Smagorinsky-Lilly模型相同，只是模型常数有区别。在流动场的低雷诺数区域，上面的函数就小于零，从而只有分子粘性起作用。所以，基于重整化群理论的亚网格模型对流动转捩和近壁流动问题有较好模拟效果。

3.3.6.3 大涡模拟的边界条件

对于给定进口速度边界条件，速度等于各个方向分量与随机脉动量的和，即:

ui=+Iψu

其中，I是脉动强度，ψ是高斯随机数，满足ψ=0，=1。

如果网格足够密并可以求解层流底层的流动的话，壁面切应力采用线性应力应变关系，即:

ρuyu=η 3,67 uημ

如果网格不够细，则假定与壁面邻近网格质心落在边界层对数区内，则:

ρuyu1=lnE 3,68 uηkμ

其中，k=0.418，E=9.793。

介绍湍流 54_湍流动能

表3,1，雷诺平均模型的比较模型名字优点

Spalart-Allmaras 计算量小，对一定复杂程度的

边界层问题有较好效果标准k?ε RNG k?ε k?ε

雷诺应力模型

缺点

计算结果没有被广泛测试，缺少子模型，如考虑燃烧或浮力问题

应用多，计算量合适，有较多数据对于流向有曲率变化，较强压力梯度积累和相当精度有旋问题等复杂流动模拟效果欠缺能模拟射流撞击，分离流，二次流，受到涡旋粘性各向同性假设限制旋流等中等复杂流动

和RNG模型差不多，还可以模拟圆受到涡旋粘性各向同性假

设限制口射流问题

考虑的物理机理更仔细，包括了湍流 CPU时间长，动量和各向异性影响湍流量高度耦合。

第四节，湍流模型算例及其设置

湍流模型设置命令:Define-model-viscous

???

无粘，层流和湍流 ?

???????

湍流模型选项 ?

???????

近壁处理方法选择 ?

?????

附加湍流选项 ?

算例分析:有换热的腔道流动问题

adiabatic wall

cold air T = 0 ?constant temperature wall T = 100 ?F

10 ft

算例二，圆柱绕流

Compute drag coefficient of the cylinder

步骤:

1，确定雷诺数，ReD=24600

2，钝体绕流，后面有不稳定的涡旋脱落。采用RNG k?ε模型，壁面处理是双层区模

型; 3，网格处理:近壁网格加密，由于是双层区模型，需要网格划分到y+=1

范文四：关于边界层湍流能量耗散率的研究

关于边界层湍流能量耗散率的研究

第29卷第1期

2011年1月

河南科学

HENANSCIENCE

V01.29No.1

Jan.2011

文章编号:1004-3918(2011)01—0006—04

关于边界层湍流能量耗散率的研究

郝芹

(华北水利水电学院水利职业学院,郑州I450008)

摘要:应用PIV系统对湍流边界层内的耗散率标度律进行了实验研究,结果显示,对不同尺度上和3个不同法向

位置湍流耗散率标度律的湍流耗散主要发生在小尺度上,即湍流耗散率标度律在小尺度上具有普适性.并且随着尺

度的不同,近壁面小尺度上的湍流耗散比较明显,远壁面上的湍流耗散率较小.

关键词:标度律;能量耗散率;PIV

中图分类号:0357文献标识码:A

湍流中多尺度结构包含着关于湍流运动非常丰富的信息,要对湍流

中的脉动机理进行深入研究,就需要

将湍流脉动信号分解成不同尺度的运动进行深入研究.这对于正确深入理解湍流产生,维持,发展和演化

的物理机理,发展湍流理论,建立符合湍流物理机理的工程湍流模式,解决重大工程技术(如航空航天工程和

水利工程)中的湍流问题具有重要意义和应用价值.

根据文献…平板湍流边界层近壁区域的平均能量耗散率可由下式得出121:

(v+):A堕:堕一.(1)

Yz,+J

这里=l=2.5,通过量纲分析可以将Kolmogorov耗散尺度与平均能量耗散率联系起来:

+)-()l,4:cy.(2)

根据标度律的公式,湍流能量耗散率的P阶结构函数[3】为

(r)=<>.C.(3)

则湍流能量耗散率的标度指数为:

lnr’(4)

这里的湍流耗散率是利用轴对称假设得出的,即

[8((+2(()2)+2((警2)一((j2)],

是湍流的能量耗散率:

fdU\

rari+l_-_J’

湍流耗散率的传统测量方法足用热线风速仪(HWA)测量出相关数据,进而计算出湍流的耗散率,根据

其实现的方法同可以通过以F三种途径:从时间序列获得能谱,空间关联系数及商接测量脉动速度的空间

梯度f等二种方法获得湍流的耗散率.目前大多学者是应用探针直接测量脉动速度的空间梯度,通过泰勒

假设,及湍流公式来得到湍流耗散率.Tsinobertj运用专门设计的具有21根热线组合探针,通过其定义的湍

流耗散率公式得到了耗散率.Andreopoulos与HonkanF]采用了9根热线的组合探针通过简化假设得到了

收稿日期:2010-07—12

作者简介:郝芹(1975一),女,河南郑州人,讲师,硕士研究生,主要从事工程结构的非线性动力学研究.

郝芹:关于边界层湍流能量耗散率的研究7

充分发展管流内耗散率分布.这种方法虽然可以对耗散率定义式中各项进行直接测量,但是由于多热线探

针自身的复杂性,使探针空间分辨率受到了限制.同时随着探针结构的复杂化,其本身对流场的干扰以及

各热丝间的干扰也将给测量带来更大的系统误差阎.本文应用PIV测量技术,相对于别的测量方法其最大

的优点是实现了同时多点测量的”平面测量”.

1实验设备及相关参数

本文主要是使用粒子速度测速仪(PIV)来对实验段进行测量,具体的实验设备由PIV系统和实验水槽

组成.PIV系统由CCD相机,激光器,同步器,帧采集器以及计算机构成.

为降低来流湍流度,驱动水泵位于实验段下游,在实验段上游收缩段前设置了5cm长的方格蜂窝器,实

验段上部开口,形成自由表面,以满足零压梯度要求.

实验是在y平面内进行的,动量损失厚度雷诺数为Reo=2694,速度矢量场大小为92x68,即方向92个

矢量,Y方向68个矢量,质询窗大小为16xl6像素;对应流场的实际大小为l0.7mmx7.9mill,采集速度为每

秒7.5幅矢量图,速度计算仍然采用FFvr相关方式进行.矢量图总数为2500幅.

2实验数据分析

为了对湍流标度指数进行详细的研究,这里对不同的尺度a=2dx,a=8dx,o=l6dx和a=32dx处耗散率标

度指数随阶数P的变化情况进行研究,同时也对不同阶数p=l,p=2,p=3和p=6时不同法向位置y+=9,y+=67

和y+=113的耗散率随尺度变化进行研究.

?=2dx

l2345678

n=16dx

a=8dx

图1耗散率在不同尺度不同法向位置处标度指数随p的变化情况

Fig.1TilerelationexponentsofdissipationlateandPascomputedindifferents

calesandwalldistances

4321O?23oooo加加

765432O

OOOOOOO

一

8河南科学第!’)第1期

由图1可以看出,在a=2dx时,耗散率标度指数表现很好的线性关系且直线的斜率为负值,随着P的增

大逐渐减小,随着尺度的增大,湍流耗散率标度指数随P变化的趋势发生变化,随着P的增加,耗散率标度指

数减小的幅度变小,当a=16dx时,呈现先增加后减小的趋势,当a=32dx时,呈现增加的趋势.且在y+=113

的位置表现的尤为明显.这也说明了湍流能量的输运过程,大尺度湍流输出能量,而小尺度湍流则耗散能

量,随着尺度的增大,湍流耗散率标度律表现出非普适性.湍流耗散率标度律在小尺度上具有普适性.

图2给出了不同阶数耗散率ln.,(r))随尺度ln(r)的变化情况.由图可以

看出,对于一阶耗散率,随着

尺度的增加,不同位置的耗散率变化不大,基本呈一直线,也就是说它的斜率都基本为零;随着阶数的增大,

湍流能量耗散率随尺度的增大而逐渐减小,这也说明了湍流能量耗散率在不同尺度上的数值变化较大,随着

阶数增加它们的变化趋势逐渐显现出来,其中在十=113位置的能量耗散率减小得比较快,而在y+=9位置的

能量耗散率减小得最慢,且尺度越大减小得越快.这也说明距离壁面越近,尺度越小,湍流能量耗散越明显,

距离壁面越远的地方,湍流耗散率越小,且随着尺度的增加而迅速减小.这也和在近壁区湍流的间歇事件

较多,相干结构较活跃是一?致的.

一

0.5

—

一

1.5

一

三一2

2.5

—

图2不同阶数耗散率结构函数在3个不同法向位置处随,的变化情况

Fig.2Therelationstructurefunctionofdissipationrateandscalerarecompute

dindifferentorderr

3结论

本文主要对湍流能量耗散率标度律进行了研究,结果显示,随着壁面距离的增大,平均耗散率逐渐减小,

这和近壁面湍流的相干结构较为活跃,距离壁面越远湍流越接近各向同性是一致的;而对不同法向位置不同

阶数的耗散率随尺度的变化作了研究,结果发现近壁面小尺度上的湍流耗散比较明显,远壁面上的湍流耗散

率较小.

2011年1月郝芹:关于边界层湍流能量耗散率的研宄一9一

参考

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

文献:

LandauL,LifshitzE.MecaniquedesFluids[M].Mir,Moscow,1973.

杨宇.湍流边界层多尺度间歇性检测和控制[D].天津:天津大学,2004.

KolmogorovAN.ThelocalstructureofturbulenceinincompressibleviscousfluidforverylargeReynoldsnumber[J].Dolk

AkadSSSR,l991,434:15一l9.

赵伟.湍流边界层中相干结构及统计量sL标度律的研究[D].武汉:华

中科技大学,2006.

EisnerJW,ElsnerW.OnthemeasurementoftheturbulenceenergydissipationLI].measSciandTech,1996,7:1334—1348.

TsinoberA,KitE,DracosT.Experimentalinvestigationofthefieldofvelocitygradientsinturbulentflows[J].JFluidMech,

1992,242:l69—192.

AndreopoulosY,HonkanA.Experimentaltechniquesforhighlyresolvedmeasurementsofrotation,strain,anddissipationrate

tensorsinturbulentflows[J].MeasurementScienceandTechnology,1996,7:1462-1476.

王汉封,柳朝晖B福水,等.用PIV数据估算槽道内湍流动能耗散率[J].

化工,2004,55(7):1066—1077.

StudiesonEnergyDissipationRateinTurbulentBoundaryLayer

HaoQin

(NoahChinaUniversityofWaterConservancyAndElectricPower,Zhengzhou450008,China)

Abstract:InthethesisPIVsystemisusedtomeasurethevectorinordertoresearchenergydissipation.Experimental

resultExperimentalinvestigationindicatesthatthescalinglawofenergydissipationrateisuniversalatsmall

scale.Theturbulentdissipationismoreobviousinsmallscalethaninlargescalewithdifferentwalldistance.

Keywords:scalinglaw;energydissipationrate:PIV

范文五：利用风廓线雷达谱宽反演晴空湍流耗散率

() 文章编号 : 100124322 20081021608207

利用风廓线雷达谱宽反演晴空湍流耗散率

1 , 2 1 1 涂爱琴, 董德保, 翁宁泉

( 1 . 中国科学院安徽光学精密机械研究所大气光学中心 , 合肥 230031 ; 2 . 山东省气象局大气探测技术保障中心 , 济南 250031)

摘要 : 从晴天风廓线雷达的回波信号机制出发 ,研究了利用风廓线雷达谱宽反演晴空湍流耗散率的方

法 ,即从雷达测量的多普勒速度谱宽中提取出湍流因素引起的谱宽 ,进而由湍流谱宽和湍流动能耗散率的关系

计算出湍流动能耗散率。通过对合肥地区风廓线雷达谱宽资料的分析 ,得到合肥上空 300 m 高度上湍流动能

- 4 - 1 2 3 耗散率的大小分布在 10 ～10 m/ s之间 ,不理论相符 ,说明利用雷达谱宽反演湍流耗散率的方法可行。

关键词 : 大气光学湍流 ; 耗散率 ; 谱宽 ; 湍流动能 ; 风廓线雷达

中图分类号 : P412 . 25 ; P425 . 2 文献标志码 : A

光波在通过大气传输时会因产生抖动、闪烁等一系列湍流效应而破坏光波参数 ,因此在研究光波和电磁波

[ 123 ] [ 4 ] 在大气中传输时 ,大气光学湍流是非常重要的。风廓线雷达由于其回波中含有丰富的湍流信息,且具有

[ 526 ] 时空分辨力高和可以无人值守连续工作等优点,成为湍流研究最有潜力的测量工具。本文介绍了使用谱宽

法反演湍流动能耗散率的方法 ,并挑选 Ai r da16000 风廓线雷达在合肥地区探测的晴空谱宽数据对方法进行了

验证 ,在此基础上分析了合肥地区湍流动能耗散率的特征。

1 仪器介绍

() 研究所用设备是北京爱尔达电子设备有限公司生产的 Ai r da16000 对流层风廓线雷达简称 Ai r da16000。该设备自 2006 年 7 月 15 日在中国科学院安徽光学精密机械研究所安装调试成功以来 ,一直连续工作 ,迄今已

[ 7 ] 经积累了近两年的宋整雷达数据。表 1 列出了 Ai r da16000 的主要技术参数。

表 1 Airda16000 型对流层风廓线仪的主要技术参数

Table 1 Ma in para meters of the Airda16000 tropospheric wind prof iler

t ra nsmit ted number pea k p ul se widt h p erio d of range gate spacing bea m mo de ()widt h/ ? μμf requency/ M Hz of bea ms po wer/ kW of lo w mo de/ m of lo w mo de/slo w mo de/s

co mpo und mo de , 450 5 3 25 . 6 0 . 8 65 . 3 150 high/ lo w

number of number of F F T ra nge gate number of number of p ul se widt h of p erio d of ra nge gate s of point s of lo w spacing of range gate s of F F T point s of μμ high mo de/s high mo de/s lo w mo de mo de high mo de/ m high mo de hi gh mo de

18 256 8 . 6 120 . 3 800 22 256 研究方法和结果分析2

2 . 1 提取湍流引起的谱宽

2 . 1 . 1 影响速度谱宽的因子

多普勒速度谱宽是由雷达体积内速度场的时间和空间变化引起的 ,即速度谱宽描述的是有效照射体积内

[ 8 ] ( ) 散射粒子相对运动的剧烈程度。一般来说 ,影响速度谱宽的因子包括以下几个方面: 1大气的湍流运动 ; () () () () () 2风切变 ; 3波束宽度 ; 4数据处理 ; 5镜式反射 ; 6不同直径的降水粒子产生的下落末速度的不均匀分() 布 ; 7重力波。

由于镜式反射通常只和稳定的大气层有关 ,在研究湍流时常有意忽略它对谱宽的变窄效应 ;而在晴空大气条件下 ,降水粒子的下落末速度和重力波的影响也可以忽略。因此 ,在研究晴空大气湍流时 ,只需考虑如何从

3 收稿日期 :2008204215 ; 修订日期 :2008209216

( ) 作者简介 :涂爱琴 1983 —,女 ,安庆人 ,硕士 ,助工 ,从事风廓线雷达遥感研究 ; bet r ue_t u @so hu . co m 。

[ 9 ] () ( ) 测量谱宽中消除 2～ 4的影响就可以了。根据 Do via k 和 Zr nic 的研究,影响谱宽的各影响因子彼此独

2 2 22 2 2 σσσσσσ立 ,于是可以假设测量谱宽是各影响因子引起的谱宽之和 ,即 :m =sh +b+ +,其中是雷达测m 量 t x ea s ear eam ea s

2 2() () σσ的谱宽速度方差为了计算方便 ,都用方差表示谱宽,表示风切变引起的谱宽 ,表示波束宽度引起的 sh b ear eam

2 2 σσ谱宽 ,t表示由于大气湍流引起的谱宽 ,表示由x 于数据处理引起的谱宽增宽值。如果能够分别计算出风切变 ,波束有限宽度和数据处理引起的谱宽 ,并将它仧从雷达测量的谱宽中消除 ,就能够得到大气湍流因素引起

的谱宽值。

2 . 1 . 2 风切变和波束宽度引起的谱宽

当平均风存在切变时 ,雷达脉冲体积内平均风速和风向会有很大的变化 ,从而导致多普勒速度谱宽产生一个较大的增量。此外即使风切变为零 ,由于波束存在一定的水平宽度 ,不波束轴线垂直的横向风在偏轴线方向上就有径向分量 ,也会导致多普勒速度谱展宽。因此 ,在使用谱宽分析大气湍流时 ,必须要对风切变和波束宽度引起的谱宽展宽值进行修订。

[ 10 ] 根据 G. D . Na n st ro m 的研究,风切变和波束宽度引起的谱宽的计算公式为2 2 θθ 5 u 20 . 52 2 0 . 52 2 2 α α( α) σσ uco s -si n uR co s+= + 0 0 0 shear beam3 3 5 z 22 2 θθ0 . 5 5 u 0 . 5 5 u ΔR 2 2 2 2 2 ( )( α)( α) ( ( α) αα( 1 ) )[ 3 + co s 4- 4co s 2] R0+ [ co s 4+ si n co s ] 3 12 24 5 z 5 z

θα式中 :是雷达半波束宽度 ;是波束天顶角 ; u是高度 R处的水平风速 ; 5 u/ 5 z 是高度 R 处的水平风速垂 0 . 5 0 0 0

Δ 直切变 ;R 是距离分辨率。

图 1 是 2006 年 12 月 1 日 300 m 高度上的风切变曲线图 ,实验用雷达 Ai r da16000 的波束宽度为 3?,则由 () 式 1可求出当日对应的风切变和波束宽度引起的速度谱宽 ,如图 2 所示。由图 2 可以看出 ,垂直波束的风切变和波宽引起的速度谱宽很小 ,即垂直波束的测量谱宽受风切变和波宽的污染较小 ,因此 ,在信噪比相当的情况下 ,垂直波束比倾斜波束反演的耗散率更可靠。

Fig . 1 Grap h of wi nd shea r Fig . 2 Grap h of sp ect ral widt h co nt ri but ed by wi nd shea r a nd bea m widt h 图 1 风切变曲线图图 2 风切变和波宽引起的速度谱宽曲线图 2 . 1 . 3 数据处理引起的谱宽

2 σ x 数据处理引起的谱宽值项说明了信号处理和方差估计中所有的谱宽影响因素戒误差源。误差源大致有低信噪比和窗函数两种。针对这两种误差源 ,采取了以下方法来减小误差 : 首先 ,以信噪比等于 5 dB 为门限 ,将信噪比小于 5 dB 的数据全部滤去 ,保证参加计算的数据的质量 ;其次 ,由于为了减弱谱泄露而使用的窗函数会导致测量得到的谱宽始终不实际谱宽成比例 ,因此只要增宽比例系数确定了 ,就可以从测量谱宽直接计

[ 11 ] 2 2 σσ算出加窗函数引起的谱宽。本文比例系数直接使用 S. J aco by2 Ko al y 总结的 4 %,即= 0 . 04。因数据 x m ea s 处理引起的谱宽不测量谱宽成比例 ,所以该项只会影响耗散率的大小 ,而对其变化趋势影响比较小 ,从而可以

() 在以后的对比实验中 ,根据雷达测量的耗散率和其它常规测量工具如超声风速仦测量的耗散率的吻合程度对该系数进行不断的修正。

2 . 1 . 4 湍流引起的谱宽

2 2 2 σ(σσ风切变、波束有限宽度和数据处理引起的谱宽都能通过公式戒经验计算出来 ,则由式=- +sh t ea s m ear

2 2 σ) σ+可以很容易就求出湍流因素引起的多普勒速度谱宽 ,如图 3 所示。 beam x

1610 强激光不粒子束第 20 卷

2 22 2 (σσ σσsh +b +但是 ,当湍流很弱时 ,由式=- t m ear eam ea s

2 2 σ) σ计算的可能是负值 ,此时不进行耗散率的计算 ,而 x t

直接赋予其零值 ,表示此时的湍流是很微弱的可忽略的

小尺度湍流。

2 . 2 计算湍流动能耗散率

() 湍流谱宽一般讣为是雷达分辨体积空间滤波和多

( ) 普勒时间序列持续时间时间滤波引起的 ,这两种滤波

效应非线性的综合叠加在一起影响湍流谱宽。对于各向

同性的均匀湍流 ,湍流引起的谱宽和湍流动能耗散率的 [ 12 ]Fig . 3 Grap h of sp ect ral widt h co nt ri but ed by t ur bulence 关系为图 3 湍流引起的速度谱宽曲线图 2 2/ 3 εk1 2 - 11/ 3 2 A ( )( ) σsi nck L / 2 = 1 -k{ 1 -2 t 2 π 4μ k

2 22 2 2 ( )e xp [ - ak bk]} d kd kd k( )+ k1 1 2 3 2 1 -3

2 2 2 2 式中 : A 是 3 维科尔莫戈洛夫常数 , 介于 1 . 53～1 . 68 之间 , 实验中取 1 . 6 ; k= k1 + k2 + k3 , k是沿着波束方向 1 的波数 , k2 和 k3 都是垂直于波束方向的波数 ;积分的范围是所有波数空间。三重积分式的指数项表示空间滤

波效应 ,而 si nc 函数项表示时间滤波效应。

[ 12 ] 由于三重数值积分在运算上实现起来很复杂 ,A . B . Whit e 于 1999 年提出了一种二重积分简化模型,简

2 2 ( 化模型通过把三重积分式转化到球坐标系后用 e xp -( ) ) x/ 3近似 si ncx得到

2/ 3 2 ππ / 2/ 2 ε AL 3 2 2 2 2 2 2 1/ 3 2 Γ( )Φθθ( θθθ)( )σ3 t = 122/ 3ddsi n b co s + a si n + co s <,从而得到si n="" 0="" 0="" 12="">

3 3/ 2 - 3/ 2ε σ( π) ( )t 4 = 4/ A I

θθτ(τ式中 : a 是囿波束横截面的半径 , a = Rsi n? R; b 是脉冲长度的一半 , b = h/ 2 =c/ 2 是脉冲宽度 , c 是 0 0 . 5 00 . 5

) ( ) 光速; L = vt, v是抽样时间间隔内的平均水平风速 , t是平均周期戒停留时间 , t= N T N N + t,其 t d t d d SA FF T PA 0

( 中 N 是谱平均的个数 , T 是脉冲重复周期 , N 是傅里叶变换的点数 , N 是脉冲平均的个数相干积分个 SA FF T PA

) 数, t是傅里叶变换处理时间。A . B . Whit e 通过比较二重积分简化模型和精确三重积分的计算结果 ,发现二 0

[ 12 ] () 重积分引起的误差小于 2 %,因此本文中使用二重积分简化模型即式 3来反演耗散率。

() 图 4 就是使用式 3反演耗散率的两个例子。由于受实验条件的限制 ,无法对雷达反演的耗散率进行实时对比。因此我仧通过检验雷达反演的耗散率的量级是否不理论相符来初步验证雷达反演湍流耗散率的可靠

[ 13 ] 性。吴晓庆采用 3 种方法估算湍流动能耗散率,得出 :在大气处于稳定条件下 ,3 种方法估算的湍流动能耗

2 3 散率平均值分别为 0 . 002 6 , 0 . 006 5 和 0 . 004 1 m/ s;在大气处于不稳定条件下 ,3 种方法估算的湍流动能耗

2 3 [ 6 ] ε散率平均值分别为 0 . 005 5 ,0 . 015 0 和 0 . 026 0 m/ s。戴铁丕得出结论:在对流层和平流层底部 , 值的变

0 - 7 2 3 - 4 - 1 化范围在 10～10 m/ s之间。如图 4 所示 ,利用谱宽反演的耗散率分布在 10 ～10 量级之间 ,不理论相

Fig . 4 Grap h of t ur bulent ki netic ener gy di ssip atio n rat e ret rieved by spect ral widt h

图 4 速度谱宽反演的耗散率曲线图

符 ,因此讣为利用雷达谱宽反演湍流耗散率的方法可行。

湍流耗散率的特征分析3

3 . 1 湍流耗散率的时间变化特征研究

为了研究耗散率随时间的变化规律 ,我仧对特定高度上的雷达数据进行分析。Ai r da16000 一共有 40 个距离门 ,能提供 40 个特性高度上的耗散率变化情况 ,为了研究方便 ,这里仅挑选 2 ,18 和 22 这 3 个距离门 ,即

300 ,2 700 和 5 500 m 这 3 个特性高度进行分析。

() () 图 5 a和图 5 c两幅图依次是 2006 年 12 月 3 日和 20 日 300 m 高度上耗散率随时间的变化曲线图 ;图 5

- 1 - 4 2 ( ) ( ) b,图 5 d是对应的水平风速曲线图。从图中可以看出 300 m 高度上耗散率的取值区间是 10 ～10 m/ 3 s;耗散率总体上具有和折射率结构常数类似的日变化规律 ,即在上午 10 : 00 左右开始增强 ,下午 16 : 00 左右开始减弱 ;此外 ,将耗散率随时间的变化曲线和右边的水平风速对比还能发现 ,湍流动能耗散率和风速大小具有很强的相关性 ,当风速大时耗散率也大 ,反之当风速小时耗散率也小。

( )Fig . 5 Grap h of ho rizo nt al wi nd velocit y a nd t ur bulent ki netic ener gy di ssip atio n rat e vs ti me 300 m

( )图 5 湍流耗散率和水平风速随时间变化曲线图 300 m

图 6 是 2006 年 12 月 3 日和 20 日 2 700 m 高度上湍流动能耗散率和水平风速随时间的变化曲线图。从

- 2 - 5 2 3 图中可以看出 2 700 m 高度上耗散率的取值区间大部分是在 10 ～10 m/ s,比 300 m 高度上的耗散率小 1 个量级 ;耗散率的变化趋势和风速变化趋势相近 ,不再具有明显的中午比早晨和傍晚耗散率强的特点。

图 7 是 2006 年 12 月 3 日和 20 日 5 500 m 高度上湍流动能耗散率和水平风速随时间的变化曲线图。从

- 3 - 4 2 3 图中可以看出 5 500 m 高度上耗散率的取值比较稳定 ,集中在 10 ～10 m/ s;该高度上耗散率随时间的变化不明显 ,不风速也无明显的相关性。

3 . 2 湍流耗散率的高度变化特征研究

为了研究耗散率随高度的变化规律 ,我仧对特定时间段内的雷达数据进行分析。本文中选择 6 : 00～8 : 00

这个湍流发生转换的时段。

图 8 是 2006 年 12 月 3 日和 20 日 6 :00～8 :00 湍流动能耗散率和水平风速随高度的变化曲线图。从图中

() ( ) 可以看出下列几个特点 : 1湍流动能耗散率总的来说是随高度减弱的 ; 2在 2 k m 左右耗散率的大小达到一

- 4 2 3 () 个局部最小值 ,量级为 10 m/ s; 3在 4 k m 左右耗散率有增强的趋势 ,但是不会超过地面耗散率的大小 ;这

() 种增强区一直持续 1～2 k m ; 46 k m 以上湍流很弱 ,湍流动能耗散率的大小主要不风速有关。

1612 强激光不粒子束第 20 卷

( )Grap h of ho rizo nt al wi nd velocit y a nd t ur bulent ki netic ener gy di ssip atio n rat e vs ti me 2 700 m Fig . 6

( )湍流耗散率和水平风速随时间变化曲线图 2 700 m 图 6

( )Fig . 7 Grap h of ho rizo nt al wi nd velocit y a nd t ur bulent ki netic ener gy di ssip atio n rat e vs ti me 5 500 m

( )图 7 湍流耗散率和水平风速随时间变化曲线图 5 500 m

4 结论

综上所述 , 利用风廓线雷达谱宽反演晴空湍流耗散率的过程包括两个重要的步骤 ,即提取湍流谱宽和计算耗散率。使用中国科学院安徽光学精密机械研究所 Ai r da16000 对流层风廓线雷达测量的数据 ,经过上述方

- 8 - 1 2 3 法的运算得到的湍流动能耗散率大小在 10 ～10 m/ s之间 ,不理论相符 ,说明反演的耗散率是可靠的。

( )Fig . 8 Grap h of ho rizo nt al wi nd velocit y a nd t ur bulent ki netic ener gy di ssipatio n rat e vs height 6 :00～8 : 00

( ) 图 8 湍流耗散率和水平风速随高度变化曲线图 6 :00～8 :00

在此基础上 ,本文还以 2006 年 12 月 3 日和 20 日为例 ,分析了湍流动能耗散率的时间变化特征和高度变化特征 ,结果表明 :湍流耗散率在 300 m 高度上具有中午比早晨和傍晚强的特点 ,且耗散率的大小和风速相关 ;在 2 700 m 高度上 ,耗散率变化规律不风速变化一致 ,且耗散率的大小比 300 m 处要小 1 个量级 ;在 5 500

- 3 - 4 2 3 m 高度上 ,耗散率取值比较稳定 ,集中在 10 ～10 m/ s,该高度上耗散率随时间的变化不明显 ,不风速也无明显的相关性 ;湍流耗散率整体上随着高度的增加而减小 ,在 2 k m 左右有一个局部极小值 ,4 k m 左右开始有一个 1～2 k m 厚的增强层。

参考文献 :

[ 1 ] ( 苏毅 ,万敏. 高能激光系统[ M ] . 北京 :国防工业出版社 , 2004 :127 . Su Y , Wan M . Hi gh ener gy la ser syst e m. Beiji ng : Natio nal Def ence In2

)dust r y Pre ss , 2004 :127

( 宊正方. 应用大气光学基础[ M ] . 北京 :气象出版社 , 1990 :772110 . So ng Z F . Fo undatio n of applied at mo sp heric op tics. Beiji ng : Chi na Me2 [ 2 ]

)t eo rolo gical Pre ss , 1990 :772110

( ) ( 翁宁泉 ,曾宗泳 ,肖黎明 ,等. 大气光学湍流测量中平均时间和原始数据的筛选[J ] . 强激光不粒子束 , 2004 , 16 9: 110121105 . Weng N Q , [ 3 ]

Zeng Z Y , Xiao L M , et al . St udy of average ti me a nd so urce dat a diff erentiatio n i n mea suri ng at mo sp heric op tical t ur bulence . H i g h Pow e r

( ) )L ase r an d Pa rt icle B ea ms , 2004 , 16 9:110221105

[ 4 ] ( 塔塔尔斯基 B H . 湍流大气中波的传播[ M ] . 北京 :科学出版社 , 1978 : 25252 . Tat a r ski B H . Wave p rop agatio n i n a t ur bulent medi um. Bei2

)ji ng : Science Press , 1978 :25252

[ 5 ] Go ssa r d E E , Wolf e D E , Mo ra n K P , et al . Mea sure ment of clea r2ai r gradient s a nd t ur bulence p rop er tie s wit h radar wi nd p rofiler s[J ] . A 2

me rican M eteorol o g ical S ociet y , 1998 , 15 :3212342 .

2 ε( ) ( [ 6 ] 戴铁丕 ,潘闻天 ,朱素芳 ,等. 用常规探空资料估算 C和垂直变化初探[J ] . 南京气象学院学报 , 1997 , 20 2:2512258 . Dai T P , Pa n W T , n

2 ε( ) )Zhu S F , et al . So undi ng e sti matio n of vertical cha nge i n Cn a nd . J ou rnal o f N an j i n g I nst i t ute o f M eteorol o g y , 1997 , 20 2: 2512258

([ 7 ] 北京爱尔达电子设备有限公司. Ai r da16000 型对流层风廓线仦技术说明书 [ Z ] . 北京 : 北京爱尔达电子设备有限公司 , 2005 . Beiji ng Ai r da

Elect ro nic Equip ment L t d. Technical sp ecificatio n of Ai r da16000 t ropo sp heric wi nd p rofiler . Beiji ng : Beiji ng Ai r da Elect ro nic Equip ment

)L t d . , 2005

( 张培昌 ,杜秉玉 ,戴铁丕. 雷达气象学[ M ] . 北京 :气象出版社 , 1992 :25228 . Zha ng P C , Du B Y , Dai T P . Rada r met eo rolo gy . Beiji ng : Chi2 [ 8 ]

)na Met eo rolo gical Press , 1992 :25228

Do via k R J , Zr nic D S. Doppler radar a nd weat her o b ser vatio n s[ M ] . A merica : Acade mic Press , 1984 :458 . [ 9 ]

[ 10 ] Na st ro m G D . Doppler rada r sp ect ral wi dt h broadeni ng due to bea mwi dt h a nd wi nd shea r [J ] . A nnales Geo p h y sicae , 1997 , 15 :7862796 .

1614 强激光不粒子束第 20 卷

[ 11 ] J aco by2Koal y S , Ca mpi st ro n B , Ber na r d S , et al . Tur bulent di ssip atio n rat e i n t he bo undar y layer via U H F wi nd p rofiler sp ect ral widt h

mea sure ment s[J ] . B oun d a r y2l a y e r M et eorol o g y , 2002 , 103 :3612389 . [ 12 ] Whit e A B , L at aiti s R J , L awrence R S. Space and ti me filt eri ng of re mo t el y sen sed velocit y t ur bulence [J ] . A me ri can M eteorol o g ical S oci2

et y , 1999 , 16 :196721972 .

( ) ( [ 13 ] 吴晓庆 ,聂群 ,方强. 近地面大气湍流平均动能耗散率测量不分析[J ] . 力学学报 , 2007 , 39 6:7212726 . Wu X Q , Nie Q , Fa ng Q . Mea s2

( ) )ure and a nal ysi s of mea n at mo sp heric t ur bulent ki netic di ssip atio n rat e i n t he bo unda r y layer . J ou rnal o f M ec h ani cs , 2007 , 39 6:7212726

Retrieval of clear2a ir turbulent dissipat ion rate using spectral

width mea sured by wind prof iler

1 ,2 1 1TU Ai2qi n, DON G De2bao , W EN G Ning2quan

( 1 . L aborat o r y o f A t m os p he ric O p t ics , A n h ui I ns t i t ute o f O p t ics an d Fi ne M ec h a nics ,

Chi nese A c a dem y o f S ciences , P. O. B o x 1125 , H e f ei 230031 , Chi na ;

)2 . Ens u ri n g Ce nte r o f A t m os p he ric S ou n d i n g Tec hnol o g y , S h an don g M eteo rol o g ic al B u re a u , J i na n 250031 , Chi n a

Abstract : Fro m t he fo r matio n mecha ni sm of t he wind p rofiler rada r’s echoe s , t hi s p ap er st udied a ret rieval met ho d of t he t ur bulent kinetic ener gy di ssipatio n rate f ro m sp ect ral widt h mea sured by a wind p rofiler . Two key step s of t he met ho d are ext rac2 ting spect ral widt h co nt ributed by t ur bulence a nd calculating t ur bulent di ssip atio n rate ba sed o n t he relatio n bet ween t ur bulent di s2 sip atio n rate and t ur bulent sp ect ral widt h . The p aper al so analyzed t he data of spect ral widt h mea sured by wind p rofiler in Hef ei

- 4 a nd t he re sult i s t hat t he t ur bulent kinetic ener gy di ssipatio n rate o ver Hef ei at t he height of 300 m di st ribute s bet ween 10 a nd - 1 2 3 10 m/ s. The re sult acco rds wit h t heo r y , so t he ret rieval met ho d i s via ble .

Key words : Op tical t ur bulence ; Di ssipatio n rate ; Spect ral widt h ; Tur bulent ki netic ener gy ; Wind p rofiler rada r

转载请注明出处范文大全网 » 用PIV数据估算槽道内湍流动