90°弯管内流动的理论模型及流动特性的数值研究110精

1、第21卷第3期2004年6月 计算力学学报Chin eseJo urn alofC ompu tati on alMecha nicsVol.21,No.3June2004文章编号：100724708(2004)032031420890 °弯管内流动的理论模型及流动特性的数值研究丁珏，翁培奋3(上海大学应用数学和力学研究所，上海200072)摘 要:从三维不可压缩雷诺时均 Navier2Stokes方程出发，对90°弯曲管道内湍流 流动进行数值模拟。网格划分采用六面体网格，湍流模型为RNGkE模型,在近壁 区采用两层壁面模型进行修正，流场的计算结果与实验数据吻合较好。在此基

2、础上 本文数值研究了来流方向对流场结构和流动特性的影响。得出在弯管流场中发生 了分离现象，且随着来流侧滑角的增大 分离区范围增大。此外，随着来流从同一侧 滑角变换至同一攻角时，横截面的二次流图像中也从具有两个对称主涡变成只具有 一个主涡的现象。关键词：弯曲管道;Navier2Stokes方程;湍流模型;流动特性中图分 类号:0357.5文献标识码:A1引言弯曲管道广泛地应用于工业、农业等机械设备上，如压缩机、泵及各种类型的热交 换器,甚至核动力的管道系统中。这些形形色色的弯管实现了流体输运和热量、质 量交换等功能。一般情况下，影响流体流动性质的因素有很多，诸如弯管的弯曲程 度，流体的来流马赫数

3、Ma,流体运动方向等。在纵多因素的影响下，弯曲管道内的流 场呈现出十分复杂的流动特性。如一定来流条件下，在管壁附近形成分离区，管道 横截面上产生二次流动，这些现象不仅造成流体总压和能量的损失，而且形成的局 部障碍区域也使流动系统的阻力增大，降低了热量、质量的交换效率。因此，弯管 内的流动一直受到内流研究者的关注123较大成绩。史峰4等人利用标准k2E湍流模型研究了 Rc弯道内的湍流流动。计 算发D=1.69的90°现利用k2E模型得到的流场时均速度分布与实验数据相比，在 弯段H =45截面前，二者符合较好，即计算误差较小；当弯段H >45°存在较大的误 差。文献采用湍

4、流大涡模型，对90。弯管进行了数值研究，拓宽了湍流大涡模型在 具有强曲率弯曲管道内流体流动问题的应用。基于以往研究的基础，本文引进Yahhot和Orszag7新近应用重整化群方法提出的 RNGkE湍流模型，建立90°大曲率弯管内部流场的Navier2Stokes方程，应用交错网格系统 下的SIMPLE算法求解，并将计算结果与文献1中的实验数据进行对比。此外，针对以往研究较少 的流体运动方向对流场结构及流动特性影响，本文着重进行了分析。给出了几种来 流侧滑角和攻角下流场内重要参量的分布情况，揭示出流场中出现的物理现象和规 律。早期的实验研究以Taylor为代表，他采用激光多谱勒测速仪对

5、90°方形截面弯管内 雷诺数分别为790和40000的来流条件进行了系统的研究，给出了层流和湍流两种 流态下的时均速度分布，以及弯曲段壁面的压力分布等重要实验结果1。文献2 应用LDV测量了弯曲度为Rc D=1.69(Rc,D分别代表曲率的平均半径与弯管的水 动力直径)的矩形截面弯道内的湍流。而文献3则通过热线测速系统对三维90°弯 管内流场的气流特性进行了实验研究。在理论研究方面426,国内研究者已取得收稿日期：2002209229修改稿收到日期:20032032311 基金项目：教育部青年教师教学科研奖励基金；上海市曙光 计划项目及上海市重点学科资助项目1Navier2

6、Stokes方程时均化作者简介:丁珏(19732)，女，副研究员； 翁培奋3(19652)男，教授，博士生导师1 2数学模型 2.1控制方程 通过对连续方程和瞬间连续方程，得到直角坐标系下定常条件、不可压缩流体流动遵循的方程。=0 (i=1,2,3;j=123)xi (1)N2S方程? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All rights reserved. 第3期丁 珏，等:90弯管内流动的理论模型及流动特性的数值研究315- (© U= ® U g A ©i2uj) xjxixjxj

7、5=u,v,w,k, E其中S5,#5分别代表通用变量5的源项和有效扩散系数。J是Jacobi变换矩阵。 u,v和w分别代表曲线坐标下的速度分量。i2=1,0,i=2i 工2其中ui,p分别代表流体的平均流速和平均压力。雷诺应力张量-© uiuj使用由2.2 RNGk- E湍流模式7，也Boussi nesgue等提出的涡粘性模型来计算。RNGkE湍流模型既适应于高雷诺数情况 3数值方法3.1 离散格式,并采用SIMPLE算法来实施计算。在计采用有限体积法对控制方程进行数值离散算域内使用交错网格8。离散后方程的一般形式为Ap5p=适应于低雷诺数下的湍流流动，即提供了一个微分形式的有效

8、粘性系数表示式，来 说明低雷诺数流动效应。湍动能方程 © ujEAnb5nb+S5p其中5表示某个待求物理量,s5p为5参量方程的源项，下标nb代表控制容积P点 周围的各网格点。3.2边界条件()=EA +Gk+Gb0 k A eff5xj5xj5xj进口条件:流体均匀地进入管道中，初始运动速度为1.0m S(雷诺数40000),温度为293K。出口条件:压力出口条件。壁面条件：固体壁面上采用无滑移条件；温度和压力分别采用绝热壁和等压梯度条 件。湍能耗散率方程()O uj= A + EA eff xjxjxjC1Ek(Gk+C3E G»3C2EO2kRNGkE湍流模式只适用

9、于离开固体壁面边界一定距离的区域，因此对近壁区的湍流采用两层壁面模型的处理技巧，即将求解 区域分成:壁面流动区和核心湍流区。壁面流动区主要包括粘性子层和部分的完全 发展区。利用湍流Reyn olds数Rey=(y表示由计算域某一点到壁面的法向距A离)来区分这两个区域。当Rey>200时，采用上述的湍流模式中的k和E方程来求 解;当 Rey<200时,使用Wolfshtein 方程模型求解近壁区的湍流流动，即仍使用上 述的动量方程和湍动能k方程,但是湍流粘性系数为 A t= OC寸湍流的影响l A =cly-l E =cly1kl A 湍流的耗散其中A。k, A盼别代表k, E的反向

10、有效Prandtl数对于低雷诺数、近壁面流动情况，可以通过消除变量，导出湍流粘性系数的微分方 程：()=1.72dEA3A M-1+CMdA MA其中 M=,C M =100 A3, r =Sk E ,31 + Br-() S=(2SijSij)2, Sij=+2xjxir 0=4.38, B =0.012,C1 E =1.42,C2 E =1.683C2E=C2E +率可用E =k1.5 l计算。长度尺度I A ,l反映了壁面AA为了适应大曲率复杂弯管内对计算网格生成的要求，将上述直角坐标系下的方程变换到任意曲)，通用的输运方程为线坐标系下(N, r ,+=(5)(5)#+#+5 N J5

11、N 5 r J5 rexp(-ex p(-),)AE其中cl=?CA-0.75,A A =70,A E =2cl4算例及结果分析为了方便分析问题，将弯管分成上游直线段、弯曲段和下游直线段三个部分。 90。 弯管的方形截面的)(#5+JS5J? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All rights reserved.316计算力学学报第21卷 边长为0.04m。弯管半径比Rc D=2.3,Rc为曲率的平均半径，Rc=0.5(ri+ro),D 为水动力直径，本文中取0.04m。弯曲段内侧壁面的曲率半径ri=0.072m,

12、外侧壁面 的曲率半径r0=0.112m。上下则指均匀来流速度和它在 x2z平面投影的夹角。首 先，我们利用上述的模型和数值方法来计算0°侧滑角和0°攻角下的流场中重要物理参量的分布，并与文献1中的实验数据比较，来验证模型和算 法的正确性。图2为计算所得的在不同横截面和不同位置处,沿主流方向的时均速度分布曲线及其与实验结果的对比。可以看到，两者符合较好IS I予.1 Thr fnmAylf bmdlLtilh UtHIVl00 ti -0.4 Mk.DgI人”ffl*Jd'.c/:?i-：V14IflUflJ 4l44LAa4.1(11<HfT I r*HS.:

13、t国弋JI U I3D牛I們亠Q6ImiL"4PL1- IuriMr QJDJQA 帕Uftrq i .n# I -» I. JuuDJ' atQA140Uttjn#-ft23r了wiv"(rQP flj 04OA01lAn *<4L 环，DO04at1* WKJ计如主t(方胡绡时迪Q沁UtK垢畫网WftWttF.3dTEnArinvn uinfPvnneAri dM0.3m。弯管的几何形状如图1所示。坐标系的原点0位于。直线段I的长度均为管道入口截面的中心。流体运动方向的变化会对流场结构与流动性质产生影响。 文中,流体的侧滑角是均匀来流速度与其在x

14、-y平面投影的夹角，用B来表示;攻角A? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All rights reserved.丁 珏，等:90 °弯管内流动的理论模型及流动特性的数值研究317较好。通过上面几组数据的比较，可以说明本文建立的数学模型能够真实反映流场 中流体的运动状况，且计算结果是可靠的。4.1不同来流侧滑角下流场的流动特性现考察来流在0。攻角,侧滑角分别为0。,40 ：80条件下，管道内三维流场出现的分离 现象、横截面上的二次流图像及沿主流方向的壁面压力分布情况。4.1.1二次流图像0 °侧滑

15、角、0 °攻角来流表示均匀流体垂直进入图中xH表示沿直线段的轴向距离与水动力学直径D之比Q取0.04m。Vc为体积平均速度。r表示无量纲的径向位置，rr0)。z333管道。图4(a)为该条件下弯管的主流入口截面H=0。的速度矢量图。由于此位置受管道曲率的影响很三(r(ri-r0)为无量纲的展向位置，z3 = z z1 2;z1是小，因此没有出现明显的回流现象，二次流速度也很小。图4(b)为弯曲段出口截面 H =90°的速度矢量图。此时截面上二次流速度较大， 已出现一对涡流，且两个主涡的涡心靠近内侧壁面，说明外侧壁面附近的压力大于 内侧壁面，推挤流体向内侧壁面流动。随着流体向

16、弯管下游流去，压力梯度逐渐减小，因此出口截面位置L=0.776m(L表示沿弯管壁面中心线的沿程长度)的二次流图 像中一对主涡的涡心远离内侧壁面，接近外侧壁面，如图4(c)所示。此外，在截面的4 个顶角处,还出现了速度值较小、强度较弱的小涡，该涡的放大图如图4(d)所示TJG t 1 L - 1 r' huN f1 f 1 HH 1 I T r F 皿 Jin 1 9J -T -Ain' R it B 于 b g 3dfc H .bp p r 巾 巾 庄 w R ?：B- 5=3 r f 1 fl 1 H 1 I r n誉wr I:, FlnunVI屮、f上i' ./,r

17、E-二-FLA-干 - -riB .' - J - ' -'-厂fr 莎三二二三口M一 一讥m 二二 八”2王三二二二瓷壬；显 立三亠二二-二？右4 r殳三二二二二？ 2嗨S 二三二 mg“ 金壬=育也 r ?护亠一 二二fe匸拓匚 3-巴迁F -.- r £:二丄氏內叭託筍厂”勰总矿I湮声直密亡2 I)-.L r讥 U 勺 IJT 严；hl五 jdHuhx广二 £0弯管宽度的1 2。图3显示了计算所得的沿弯曲段外侧壁面中心线位置的压力系数cp随极角H的变化关系。cp=2(p-pref)( © Vc),pr是 H =0° ,r3=

18、0,z3=02位置处的参考压力值。观察图3可以看到,由于流体运动受到弯管曲率和离心力的影响，沿弯曲段外侧壁面的压力值呈现出先增大后 减小的变化趋势。此外，图中计算的cpH曲线与实验数据1一致性? 1995-2006 Tsin ghua Tongfang Op tical Disc Co., Ltd. All rights reserved.318计算力学学报第21图5展示了侧滑角为40。情况下沿主流方向上各横截面的速度分布情况。由 于受到来流方向和流场产生分离现象的共同影响，管道上游直线段L=0.06m位置 处，横截面上的二次流运动速度已产增加”极易导致流场中小涡消失。然而，本文计算得出的小涡

19、清晰可见，这也说明 了本文采用的数值方法是正确的。4.1.2分离现象0。侧滑角下，沿主流方向的流场分离现象不明生。图5(a)显示了在靠近后侧壁面与外侧壁面之间拐角处出现的回流区的放大 图。弯管出口截面的速度矢量如图 5(b)所示。可以看到，图中出现与图4相似的二 次流图像，即存在一对大小相近的对称1 5左右。而显。当侧滑角增大到40。时，情况已发生变化。在靠近管道入口附近处，分离现象发 生，流场中出现一个较大的回流区域，如图7(a)所示，这主要是受到流体粘性和来流 方向的影响所致。此条件下回流区的面积占弯管上游直线段部分的 且，在来流的迎风壁面侧，流体速度较大,而背风侧附近即分离现象发生的区域

20、，流体 运动速度最小。此外，在背风侧(外侧壁面)的附面层内，可清晰看到速度分布的特 点。图7(b)是侧滑角为80°时纵剖面(z3=0)上的流体速度分布图。可以看出，分离 现主涡。而且,在该截面的4个拐角处，同样也出现了强度较弱的小涡。当侧滑角进一步增大到80°时，如图6所示，由于二次流造成流体间动量和能量的交 换，使得截面上速度分布向外侧壁面和前后侧壁面推移，因此在距离外侧壁面较近 的区域又生成一对反向的次涡，如图6(c)所示。可见，此条件下管道出口截面呈现了 多个大小不一的回流区。通过以上几组图，可以得出：由于受到弯管曲率及流场逆压梯度的影响，管道截面上 产生了二次旋流。

21、该现象的出现给流场带来了横截面方向上的动能和能量交换，既改变了主流方向的平均速度，又减少了流场的总压。此外，随着侧滑角的增大，管道 出口截面出现了一对主涡、一对次涡，以及两对小涡的图像。主涡流的强度较大， 在管道存在曲率的情况下就会出现；次涡、小涡的强度依次减弱。而且，随着侧滑 角从0°增加到80°,管道出口截面上主涡流的运动速度逐渐减小。另一方面，由于数值计算中会出现寄生的 数值湍流”而使湍流粘性象十分严重。即相比于40。侧滑角情况，该条件下回流区域的范围明显增大。在发生分离现象的流场中，逆着主流方向流动的流体在来流介质的冲击、带动 下，又重新进入回流区。这一现象的发生，

22、会使流场产生局部流动障碍，增加了输运 系统的阻力。同时，使流动总能量部分以热量形式耗散掉，降低了管道传输流体的 效率。如果涉及的是热交换器的弯管，就有可能造成管内外热量交换效率降低的后 果。因此，针对本文所讨论的情况，为了实施流动控制，应尽量减小来流速度的攻角， 即均匀来流垂直弯管的入口截面进入管道时片IMH'l'讪*帖- '11|血lirSciuu 生 F-T * 隹冷4 :： W 1?料代鬻北* '、 W :，:、 1T$甲克乜A A c i-二f-i M Ai.咼P忻占嘴4 ijjv, t匸< -严11 .kTiTP-p I ' iP i h

23、ll |- y? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All rights reserved.y”f- - /'W1- -ll_ &ii：-I 1 *5-匸"、rV-S-h 士 二-丄_l_-L. - I二u阳戶 n'! i 71-1. ： -1i t,. *1 c汕huF： 士 ir 皆：h LMTJ i -a S- i山 I k+ £1Pill ,1 fl Iby' Ih- W F，- e hk- Ji -订 «:|-:r，弯管内流动的理论模型及流动特性的

24、数值研究319? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All rights reserved.320计算力学学报第214.1.3内外侧壁面的压力系数分布图8分别显示了 0°,40和80°侧滑角情况下，沿弯管壁面中心线长度的静压分布。 这里的参考压力等物理量均取弯管入口处来流的参量值（即L=0.0m位置）。图中可 以看到:对于0°侧滑角的来流图8（a）,弯曲段区域的内侧壁面（弯曲段曲率半径较小的壁面）所受压力值明显低于 外侧壁面（曲率半径较大的壁面）的压力值，这主要是流体介质沿弯管流动时受惯性

25、 离心力作用所致。此外，由于流场发生分离现象，加上流动过程中附面层变化引起 的损失，使得弯管出口截面内、外侧壁面上的压力小于进口处的压力值。图8（b）和图8（c）的曲线中，在管道入口附近的内、外侧壁面上均出现压力平台，这主要是流场分离现象发生的区域。4.2不同来流攻角下流场的流动特性上节已对0°侧滑角、0°攻角的情况进行了描述。下面，着重讨论0°侧滑角下,攻角 分别为40°,80条件下弯管内流体的流动特性。4.2.1二次流图像图9为来流攻角40。条件下弯管出口截面上流体运动速度的矢量图。可以看到，不同攻角引起截面上二次旋流的图像不同于来流侧滑角的情况。即

26、二次流图像在管 道曲率和逆压梯度的影响下，已不再生成对称的涡流，而是以一个强度较大的固体 涡形式存在。涡流运动方向仍是从外侧壁面流向内侧壁面。同时，在壁面顶角处均有4个流动方向与固体涡相反的次涡。_ 二 f*電当来流攻角增大到80°时，在管道入口附近（InJh J H , iF*j J ；i til I II * II *III.理卩：山：般血旳住奥A甲氓2 型n汕r帖坍Th vrlr-H Uy 5制叮严曲訂厂和呵« 别；“河广”厂I/>1' JA.JI-t-= ='II:*八/« “ Mf亠”* WI-* 4b . J / /；F 4 H

27、1 JI.I"P fi' 存T 'r-，.弋=pF/广2 Vi 4 y h 7齐1¥" m【：U B1H紅厂怕H"童注：IIW( F”，样 b I 咗 E Tile vrSt- ：'：i.嘔 J ; I |】计扎“ jt djt-；3Ci:' hitH?inunj1L=0.1m），来流介质强烈冲击迎风侧壁面（前侧壁面)，造成此时前侧壁面所受的压力值增大，流体运动速度也随之增大。而且，在迎风 侧壁面的阻挡下，流体在靠近两侧壁面的拐角区域急剧转弯，并沿着两侧壁面向背 风侧壁面(后侧壁面)流去，如图10(a)所示。此时，二次流图

28、像中靠近后侧壁面的 拐角处分别出现一个较对称的回流区。图10(b)显示了弯管出口截面的速度矢量图，仍呈现出单涡现象，但此条件下流体的二次流运动速度小于40。攻角情况。4.2.2内外侧壁面的压力系数分布图11显示了 0°侧滑角,40 °、80。攻角条件下，沿弯管曲壁面中心线的压力分布曲 线。可以看到，在弯曲段H =0。90°区域,曲率半径较大的外侧壁面所受的压力值仍 大于内侧壁面。而且，曲线中出现一段压力平台，这主要是流场发生分离现象所 致。? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All ri

29、ghts reserved.m j 40.相10*皿漳下主賣方压#«，冲（41 PK iXWwrt wmcmutkhi* duCl bntAidlnri dirtCWn第3期丁 珏，等:90。弯管内流动的理论模型及流动特性的数值研究3215结论本文利用数值模拟的方法研究了90。弯曲管道内流体的流动特性。其中，采用雷诺时均Navier2Stokes方程和RNGkE湍流模型来描述流.Exp eriChe ng,GuoR on gweime ntali nv estigatio nforse paratedflowi nthree2dime nsio nalbe ndductJ.ActaA

30、ero nauticaetAstro nauticaS in ica,1994,15(9):109521099.(i nChi nese)4史 峰，徐 忠，马材芬.大曲率弯道内湍流数值计算体的湍流流动，同时利用有限体积法和SIMPLE算法实施计算。计算结果与实验数 据吻合得很好，说明建立的数学模型和数值方法是正确的。此外，文章还着重讨论 了以往研究较少的不同来流攻角和来流侧滑角下，流体介质的流动特性。如横截面 上产生的二次旋流、流场中出现的分离现象、迎风及背风壁面上的压力分布等，这些物理图像是符合实际流动的基本特性。而且，本文的研究结果也为揭示弯管内复 杂流动现象的作用机理及实施流动控制提供了

31、参考。与测量J.空气动力学学报，1990,8(4):4232429.(ShiFe ng,XuZh on g,MaCaife n.Turbule ntflowi ncurvedductwithlargecurvature nu mericalc ompu tatio nan dmeasureme ntJ.ActaAerod2yn amicaS in ica,1990,8(4):4232429.(i nChi nese) 樊洪明，何钟怡,王小华.弯曲管道内流动的大涡模拟J.水动力学研究与进展,2001,16(1):78283.(Fa nHon gmi ng,HeZh on gyi,Wa ngXia

32、ohua 丄 argeeddysimulatio no facurvedductflowJ.J ournalofHydrody namics,2001,16(1):78283.(i nChi nese) 许 承,郭荣伟.一种分离流数值求解方法J.空气参考文献(References):1 TaylorAMK P, WhitelawJH,Yia nn eskisM.Curvedductswithstro ngsec on darymotio n:动力学学报,1994,12(3):3632366.(XuCheng,.An effectivemethodfornumericalGuoR on gwei

33、simulatio nofsep aratedflowsJ.ActaAerod yn a2micaSi ni ca,1994,12(3):3632366.(i nChine se)velocitymeasureme ntsofdevelo pin glami naran dturbule ntflowJ.JournalofFluidsE ngin eer in g,1982,104:3502359.2 马材芬,史 峰，徐 忠.90弯道内紊流的试验研究J.工程热物理学报,1990,11(1):37240.(MaCaife n,ShiFe ng,XuZho ng.Develo pin gturbu

34、le ntflowin90°urvedductJ.JournalofEn gi neeri ngThermo physics,1990,11(1):37240.(i n7 YakhotV,OrszagSA.Re nor malizatio ngroupan alysisofturbule nce1.BasictheoryJ.JournalofScie ntificCom putin g,1986,1(1):39251.8 范维澄，万跃鹏.流动及燃烧的模型与计算M.合肥:中国科学技术大学出版社，1992.(FanWei2che ng,Wa nYuepen g.TheModela ndC

35、o mp utatio no fFlowa ndCombustio n M.Hefei:U nive rsityofChi nese)3许 承,郭荣伟.三维弯管内分离流场的试验研究J.Scie ncea ndTech no logyofCh inaP ress,1992.( inChin ese)航空学报,1994,15(9):109521099.(Xu(下转第329页)? 1995-2006 Tsin ghua Ton gfa ng Op tical Disc Co., Ltd. All rights reserved.第3期荣见华，等：一种基于应力的双方向结构拓扑优化算法329.TheES

36、Oalgorithmcatersfort opo logy op ti invo Ivetheremovalofmaterialmizati on byallow in gtheremovalofmaterialfromall partsofthedesig nsp ace.Witha pprop riatechequer2boardc on trols an dco ntrols on the nu mberofcavitiesformed,.con strai ntsthemethodca nrep roducetraditio nalfullystressedto pologies,a

37、ndhasbee napp lie di ntothe problemswithstaticstress,stiffness,dis placeme ntetc.lfthealgorithmwasrestrictedt otheremovalofsurface2 only material,the nasha peRece ntresearch(Q.M.Quern)has prese ntedabi2directio nalo ptimizatio np roblemissolved. evolutio narystructuralo ptimizatio n(BESO)methodwhere

38、bymatenalca nbeaddedtoa ndrem ovedoff.Buttherearemuchmoreoscillatio nstatesi nop timizi ngiterati onp rocessesofthisBES Omethod,itleadstolo ngcalculatio ntimefora nop timumsolutio n.lno rdertoim provethisBES Omethod,i nthis pap er,based on stressesa ndtheevoluti on arystructural op timizatio nm etho

39、d,a p rocedurefora nimp rovedbi2directio nalstructuralt op ology op timizatio nisgive n.ltisadevel opmentan dmodificati ono fthec onven ti on alESOa ndtheBESOmethod.Twoexa mp lesdemo n stratethatthe propo sedmethodca ndealwiththesolutio no scillatory phenomenon ,a ndobtai n more op timalstructuralt op ology