（流体力学专业论文）采用非线性kε湍流模型的污水岸边排放数值模拟.pdf

采用非线性k - e 湍流模型的污水岸边排放数值模拟 流体力学专业 研究生：何丹指导教师：杨骏六 由于城市现代化和工业发展，城市生活废水和工业污水都逐渐增多。避免 环境被污水污染成为环境保护的一个越来越重要的方面。水污染防治一方面需 要研究高效的污水处理方法，另一方面需要有效地控制污水向环境的排放。 过去人们多采用实践的方式来控制污水的排放，比如，到实际环境中测景 水质，与标准水质相比较，再反过来限制，调整污水排放量。这类方法虽然切 实可行，但是开销较大。当前，随着计算流体力学的发展，污水排放过程和污 染物的浓度分布已可以由数值模拟得到。使用数值模拟的方法，在实际现场测 量之前，对污水在受纳水体中的扩散便可以掌握其大致情况。这有助于减少实 际测量的花费，使得污水排放的控制更加方便。 作为污水排放数值模拟的典型问题，河道岸边侧向排放得到了广泛的研究。 污水岸边排放属于典型的对流传质问题，且在现实情况下极大可能是湍流问题， 故湍流模型的选取在数值模拟中起着重要的作用。在目前多数研究中，采用的 湍流模型为标准k f 湍流模型，这种模型虽已被证实是一种适应范围较火的模 型，但其并不适用于流体中雷诺正应力占重要作用的情况，目前所知，具有回 流区和二次流的一些流动便属于此种情况。为了使数值模拟能够更加适应真实 条件，本文采用了一种更优越的湍流模型，非线性k ，s 模型进行计算。同时， 采用有限体积法将控制微分方程离散为代数方程，采用s i m p l e 方法求解压力。 最后，得到了合理的数值计算结果，并与实验数据和标准k e 模型的计算结果 进行了对比。 关键词；污水排放，数值模拟，非线性k 一湍流模型 n u m e r i c a ls i m u l a t i o nt ot h ew a s t e w a t e r ss i d ed i s c h a r g e s u s i n gn o n l i n e a rk - f t u r b u l e n c em o d e l f 1 u i dm e e h a n i c s g r a d u a t es t u d e n t ：h ed a n s u p e r v i s o r ：y a n gj u n l i u t h ew a s t e w a t e rf r o mt h ec i t yr e s i d e n c e sa n dt h ei n d u s t r yf a c i l i t i e sg e t sm o r e a n dm o r et h e s ey e a rb e c a u s eo ft h ec i t ym o d e r n i z a t i o na n dt h ei n d u s t r yd e v e l o p m e n t t h u s ，p r e v e n t i n ge n v i r o n m e n tf r o mb e i n gp o l l u t e db yw a s t e w a t e rg e t st ob eam o r e a n dm o r ei m p o r t a n ta s p e c to fe n v i r o n m e n tp r o t e c t i o n i no r d e rt od ot h i s ，t h e p r o b l e mo no n eh a n d ，i st os t u d ya b o u tt h ee f f i c i e n tw a y st ot r e a tt h ew a s t e w a t e r , a n do na n o t h e r , i st oe f f e c t i v e l yr e g u l a t et h ea m o u n to fw a s t e w a t r sd i s c b a r g i n gt o c a v i r o n m e n t i nt h ep a s t ，p e o p l eu s e dt oc o n t r o lt h ed i s c h a r g eo fw a s t e w a t e ri nt h ep r a c t i c a l w a y s ，f o re x a m p l e ，m e a s u r et h ew a t e rq u a l i t yo fr e c e i v i n ge n v i r o n m e n t ，c h e c ki t w i t ht h ec r i t e r i o n ，a n dt h e nc o m eb a c kt oc o n t r o lt h ed i s c h a r g i n ga m o u n t t h i sk i n d o fm e t h o d si sp r a c t i c a l ，b u ti t sh o l d i n gf a i d yl a r g es p e n d i n g n o w a d a y s ，w i t ht h e d e v e l o p m e n to fc o m p u t a t i o n a l f l u i d d y n a m i c s ，t h ep r o c e d u r e o fw a s t e w a t e r d i s c h a r g ea n dt h ec o n c e n t r a t i o nd i s t r i b u t i o no ft h ew a s t ec o u l db ee a s i l yg o tb y n u m e r i c a ls i m u l a t i o n b yt h i sm e a n s ，ag e n e r a li d e aa b o u tt h ed i f f u s i o no f w a s t e w a t e ri nt h er e c e i v i n gw a t e r sc a l lb eg o tb e f o r et h ea c t u a lm e a s u r e m e n t ，w h i c h h e l p s t oc u td o w nt h em e a s u r e m e n ts p e n d i n ga n dm a k e st h er e g u l a t i o no f w a s t e w a t e rd i s c h a r g em o r e c o n v e n i e n t a sat y p i c a l & s ei nt h ef i e l do fn u m e r i c a ls i m u l a t i o nt ot h ew a s t e w a t e r d i s c h a r g e ，t h ep r o b l e mo fs i d ed i s c h a r g ei n t oar i v e ri sw i d e l ys t u d i e db yr e s e a r c h e r s i ti sap r o b l e mo fc o n v e c t i o nm a s st r a n s f e ra n dc o u l db em o s tp r o b a b l yap r o b l e mo f t u r b u l e n c ei nr e a lc o n d i t i o n s ot h et u r b u l e n c em o d e la c t sa ni m p o r t a n tr o l ei nt h e n u m e r i c a ls i m u l a t i o n i nm o s tc a s e si nc u r r e n tr e s e a r c h e s ，t h es t a n d a r dk s t u r b u l e n c em o d e li sa d o p t e d a l t h o u g l lt h es t a n d a r dk 一# t u r b u l e n c em o d e lh a s b e e np r o v e dt ob eam o d e la d a p t i n gal a r g er a n g eo fc o n d i t i o n s ，i td o e sn o tf i ti n w i t ht h ec a s e sw h e r et h en o r m a lr e y n o l ds t r e s so ff l u i di so fg r e a ti m p o f t a n c e a s g e n e r a l l yk n o w n ，t h e s ec a s e si n c l u d ea t l e a s ts o m ec a s e sw h e r er e c i r c u l a t i o no r s e c o n d a r yf l o wt a k e sp l a c e i no r d e rt om a k et h es i m u l a t i o na d a p t i v ei n ag r e a t e r d e g r e et 0t h er e a lc o n d i t i o n s ，i nt h i sd i s s e r t a t i o nam o r es u p e r i o rt u r b u l e n c em o d e l k n o w na sn o n l i n e a rk sm o d e li s a d o p t e df o rc a l c u l a t i o n t h e f i n i t ev o l u m e m e t h o di su s e dt ot r a n s f o r mt h ec o n t r o l l i n gd i f f e r e n t i a le q u a t i o n si n t ot h ea l g e b r a o n e s ，a n dt h es e m i i m p l i c i tm e t h o df o rp r e s s u r e - l i n k e de q u a t i o n s ( s i m p l e ) i su s e d t os o l v et h ep r e s s u r e ar e a s o n a b l en u m e r i c a ls i m u l a t i o nr e s u l ti sg o t ，a n dt h e ni s c o m p a r e dw i t ht h ee x p e r i m e n t a ld a t aa sw e l la st h e r e s u l tg o tb yu s i n gs t a n d a r d k m o d e l k e y w o r d s ：w a s t e w a t e rd i s c h a r g e ，n u m e r i c a ls i m u l a t i o n ，n o n l i n e a rk - et u r b u l e n c e m o d e l 1 工程背景及研究现状 1 1 水污染的危害及我国水污染的现状 随着国民经济的发展和工业化，城市污水和工业废水的排放量将会逐渐增 多，水污染如果得不到有效的控制，将会带来各方面的危害。 水污染的危害表现在各个方面：一，直接影响饮用水的水质，危害人体健 康。二，降低农作物的质量和产量。许多事实表明，若含有有毒有害物质的废 水污水污染了农田土壤，不仅造成农作物枯萎死亡，且对土壤造成的大面积的 污染难以治理，将造成长期的损失。三，制约工业的发展。由于许多工业需要 利用水作为原料或洗涤产品，水质的恶化直接影响产品的质量，水质恶化也会 造成冷却水循环系统的堵塞腐蚀和结构问题。四，影响渔业生产的产量和质量。 水污染会造成鱼类大面积死亡，还会使鱼类和水生物发生变异。此外，水污染 还会使有害物质在水生物体内积累，减低其食用价值。五，加速生态环境的恶 化，对环境影响十分明显，造成巨大的经济损失。据世界银行估计，我国大气 污染与水污染所造成的经济损失是每年2 4 2 2 亿美元，占我国国民生产总值的 3 5 。因此一个国家在经济发展的同时，如何保护水源水质的清洁是一个根 本性的问题，也是现代人类所面临的一个值得着力研究和解决的菇同问题。 目前，我国的水污染状况主要涉及到有机污染、重金属污染、富营养污染以 及这些污染形式共存的复合性污染。应该注意的是，受到污染的河流往往同时 存在大量悬浮物，它们中的相当一部分是有机物，排入水体后先是沉淀至河底 形成沉积物。沉积物是水体的一个潜在污染源。近年来难降解合成有机物污染 受广泛注意，它们即使在十分低的含量的情况下也可能对人体健康产生直接危 害，如致癌、致突变等，严重影响了人民群众的生活和工农业生产。我国从上 世纪7 0 年代中期开始着手进行水污染防治已经3 0 年了，但到目前为止国内一些 江河湖海的污染情况却仍然十分严重，令人担忧。根据( 2 0 0 2 年中国环境状况 公报，我国主要水系的水质仍不能达到其功能的要求，七大水系7 4 1 个重点监 测断面中，2 9 1 的断面满足一至三类水质要求，3 0 的断面属四至五类水质， 劣五类水质断面占4 3 6 。我国污水、废水排放量逐年增大，据统计，城市生 活废水约占废水排放量的4 0 ，工业废水占6 0 。而近年来由于工业的快速增 长、人口压力以及农药化肥使用量的不断增加，我国地面水、地下水的质量都 有较大的下降，全国估计每年水污染造成的经济损失约4 0 0 , ( l 元，因此保护水资 源、防治水体污染已成为我国政府十分重视的大问题。 1 2 水污染防治与污水岸边排放数值模拟的意义 对于水污染的防治，首要的工作是开发和采用适宜的城市污水和工业废水 处理技术，对不同性质的废水进行有效的处理。在这方面我国随着科技攻关 项目“水污染防治及城市污水资源化技术”的展开，在吸取国外经验，并结合 自身国情的基础上，通过发展充分利用大自然净化能力的稳定塘、土地处理、 深海排放等技术；能实现污水回用的深度净化技术：使用量大面广的造纸、印 染、高浓度有机废水治理技术，逐步搭建起了自己的废水处理体系。但是另 方面，处理过的废水排入水体，如果不控制排污量和采用适当的排入方式，污 染物不能有效扩散，也会对水环境造成污染，对排放废水的污染物浓度进行控 制可以采用实时监控的方式。目前对于一些重大污染源，国家投资建造了相应 的远程监控系统。但此类系统投资大难蛆普遍推广，对于大量城市污水和工 业废水的排放源，应该寻求更经济实用，更具普适性的污水岸边排放控制措施。 近年来，随着湍流模型网格技术的发展和计算机速度的加快，对现实中 复杂流场进行传质传热数值计算已经实用化，通过计算机进行数值模拟，对于 不同的流域边界，排放口条件下的污水岸边排放，均可快速得到浓度分布的详 细信息，配合现场实验，可以为排放的可行性分析，排放口的最优工艺参数选 取等提供合理的指导。近年来，污水岸边排放的数值模拟已逐渐得到重视和发 展。 1 3 污水岸边排放数值模拟研究现状 首先需要说明的是，污水岸边排放数值模拟一般按模拟对象可分为针对排 放口近区流场或针对排放口远区流场两种情况。关于排放口远区流场的数值模 拟，可参见文献【l l ，本文主要讨论针对排放口近区流场的数值模拟。 在近区( 或称混合区) ，污染物与环境水体剧烈掺混，形成结构复杂的回流 区，各个特征量变化强烈，故一般的数值模拟方法需要在对流传质问题的基本 方程组的基础上，配合适当的湍流模型进行流场的求解。在河道水深与水宽尺 方程组的基础上，配合适当的湍流模型进行流场的求解。在河道水深与水宽尺 度相当或者河床起伏较大的情况下，直接进行紊流模型的三维求解【2 1 。而在深 宽比较小的情况下，可以采用深度平均处理，将问题从三维形式化简到二维， 进行二维求解p i 。深度平均处理是将三维问题的微分方程在深度方向上积分， 将关于原变量的三维问题转变为关于深度平均量的二维问题的一种化简方法。 上世纪7 0 年代，国外的研究者们对侧面注水的槽道流这样的基本问题进行了研 究州，他们从这些基本研究中创立了深度平均理论，后来成为了宽浅水道岸边 排放流场计算的理论基础。深度平均的处理方法能大大化简计算，且大量实际 问题符合于作该处理的前提，是目前在实际研究中得到广泛应用的方法。近期， 污水岸边排放数值模拟研究更侧重将理论联系实际，研究的兴趣从规则形状河 道的侧面排放这样的条件比较理想的问题转移到解决更加实际的一些问题上， 比如，发展贴体坐标和网格技术以将数学模型应用范围扩展到能适应自然河道 的计算【6 1 ；发展二维动态水质模型及三维稳态水质模型通过计算非定常流来预 报污染物扩散混合的过程【7 1 ；模拟实际河道排放口附近污染物运动轨迹【8 】，等 等。 1 4 本文的工作 如上所述，从查阅到的资料看来，近期污水岸边排放数值模拟研究更注重 将数值模拟深入到具体的实际问题当中，去解决现实中迫切需要解决的各类实 际问题，而污水岸边排放问题基本的流场计算模型和计算方法已经比较固定。 虽然如此，对于具有规则几何条件的岸边排放近区计算这样的基本问题，若能 提高其精度仍是有其意义的，污水岸边排放近区流场本身流动较复杂，加之近 区流场的出口又是远区流场计算的进口，故对近区流场若能进行更精确模拟， 必将有助于寻求提高其他实际的岸边排放问题模拟精度的方法。 目前岸边排放近区流场数值模拟中较多使用的是深度平均的标准k 湍流 模型，标准k 湍流模型是一种被证实具有较高精度和较广泛适应性的模型。 但对于存在各向异性，二次流或有较大的流线曲率的复杂流场，标准k 湍流 模型的计算结果将在某些方面失真，而在天然河道中，这些情况可能是存在的。 故若能引入更能适应流场复杂性的湍流模型，则有益于提高计算的精度。目前， 已有研究者将应力通量代数模型应用到求解热水岸边排放近区流场 1 0 1 ，应力 通量代数模型为二阶封闭紊流模型的一种简化形式，在文献【1 0 中，采用此模 型得到了比标准k f 模型更真实的结果。本文的工作也在于通过采用新的湍流 模型提高污水岸边排放数值模拟的精度。考虑到计算精度和开销的平衡，本文 引入了一阶封闭模型中适应性较广且便于实施的非线性k 模型对流场进行了 计算，通过所得结果与文献 3 】中的标准k s 模型计算结果及实验结果的对比， 证实了非线性k 模型对流场计算的改进作用。另外，非线性k 模型计算耗 费的机时也仅比标准k 一模型情形略有增加，小于二阶封闭类的模型，故其可 作为污水岸边排放近区流场计算选用的一种新的湍流模型。 对流传质问题基本的数学模型是确定的。数值计算方法上，有限体积法及 速度压力修正( s i m p l e ) 的方法也被大量研究和实践证实为处理此类问题效 果较好的方法，得到普遍采用。本文将在第二，三章中简要介绍问题数学模型 的提法和相应的数值计算方法。在第四，第五章中，引入非线性k 一6 湍流模型， 对问题分别使用标准k 一湍流模型和非线性k s 湍流模型进行计算，并对结果 进行比较和分析。 对上述采用标准k s 湍流模型和非线性k s 湍流模型的污水岸边排放问题 有限体积法求解过程，本文使用c + + 语言编制了程序，并使用了v i s u a lc + + i d e 编制了w i a d o w s 程序对计算结果进行了绘图。 2 污水岸边排放问题流场计算的数学模型 2 1 控制方程组 在河道的深度与宽度同尺度的情况下，污水岸边排放问题属于三维问题， 但对于实际中存在的大量宽浅河道，可以通过水深平均处理把问题化为二维。 考虑湍流的对流传质问题基本方程组在直角坐标系下有如下形式： 螋。o a z 二 毒c p 叩小詈+ 毒h 等+ 等) + 毒( _ p 丽) 毒伽- e ) 。击卜詈一p 币) 。 u 掣一毒( 尝等) + g p s q 掣一旦a x , f k 堕a , 旦a x , 1 + c 。詈g c ：。p i f , 2 式中，u 为时均速度，标量0 为污染物的质量浓度( 时均值) 。 k ，s 为湍动能及其耗散率。 c 1 ，和c ：。为湍流模型常数。 “为脉动速度，因此一u i u ，为雷诺应力。同样，疗为浓度的脉动值，u ，口为浓度 的湍流通量。 角标i ，j 取值从1 到3 ，并遵从求和约定。 在本文中，均用大写字母u ，0 等表示物理量的时间平均量。而用小写字 母“，口等直接表示脉动量，不再加撇。这样的记法方便文中介绍各种湍流模 型时，将会进行的大量与脉动量相关的叙述。 将上式在水深方向上积分并化简可以得到各深度平均量的方程组： 5 竺竺。0 a 茸a ， 掣+ 掣一去( 弘警) + 专詈) 一p ( 警+ 警卜昭芸 缸眇缸l 缸砂i毋l 批妙j 仳 + 字+ 丢捌p 妙一u 池+ 丢泓p 妙一u 肌一咄 掣+ 掣一去( 尝) + * 必a y ) 一p ( 警卫a y 卜茜 缸砂 缸缸7 妙i 。7 i 以j d y + 曼 蔓十丢鲁p p 一u 肌一矿k + 丢舒p 帜一矿) 2 出 t a ( p u o ) + 掣一静塑a x ) + 专( 塑a yh 警+ 钟警 融 却船i 。j 妙。j l 缸妙 + 砉瓤p 移一u x o - o k + 砉泓p 帜一矿x o - e k 掣+ 掣- 跚+ 爿卦北+ 尝) 割 一p l “ - i i a u + - i va 万u + 石u ua 缸v 叶芳) _ 卵+ 暇 掣+ 掣* 水誓p tl 川a e l 十北+ 等) 爹】 一c 詈f 詈+ 面詈+ 面警+ 芳) - p 等+ 绒 ( 2 2 ) 式中，u - 知u 幽 y 一暂k 出， 6 ，三亡 d z h j o 为速度和浓度的深度平均量，u ， ，o 为未经过水深平均的原始变量，换 用这种写法是由于此后讨论的变量均为水深方向的平均量，故水深平均量直接 记为了u ，y ，0 ，而未经水深平均的原始变量改记为u ， ，而式中 上横线表示取湍流时均，不是指水深平均。( 2 2 ) 中，动量方程和浓度方程的 右边后三项，k 和s 方程的右边最后一项是由深度平均处理新产生的项。其中， k ，b 为水面剪切力，为底面剪切力，吼为通过自由面的标量通量， 这里对于浓度来说为0 。咒，只是k 和5 方程水深平均过程中垂向速度分布不 均匀造成的产生项。可按如下摩阻定律计算：f 。p ；c ，u 2 + v 2 厂， p ，c i 矿p 2 + 矿2 y ”，其中c ，为阻力系数。而p h = c 。u 3 h ， p w - c , u , 4 1 1 2 。其中，c t - 两1 ，c 。一1 3 6 c , r 2 c 厅，u ti c ，u + v 2 卜 h 为水深。水面剪切力k ，按照自由面的刚盏假设可以略去，由于浅水问题 速度及浓度在深度方向变化较小，动量方程和浓度方程后面的两项积分也可以 略去，最后我们用于求解浅水岸边排污问题紊流流场的控制微分方程组简化为 如下形式： 7 1 a ( p u _ u ) + 掣- 昙( 警) + 号卜詈) 一p ( 警+ 萼卜罢一譬 掣+ 掣= 静针* 豺p 降斟考一盟p a ( p 缸u o ! + 掣一去( 静* 豺p 降爿 掣+ 掣。北+ 尝) 讣北+ 纠剖 一p f 婴+ 茹型+ 面里+ 坐1 一胪+ 以 。i 缸 a y a x a yj ” 掣+ 掣；北+ 尝肼北+ 尜】 _ - p 詈( 万尝+ 西等+ 面豢+ v - c 面a r l j - c , 2 p 1 。- + 以 ( 2 _ 3 ) 其中，p 为压力，即流体微团3 个方向正应力的平均值，在本文的情况下等于 舟咖 这里连续性方程未再写出，因为后面使用压力速度修正计算方法时，连续 性方程便会隐含地使用到，而不需显式对其求解。 ( 2 3 ) 中的雷诺应力1 5 2 ，而，和浓度速度相关面，面需要通过湍流模 型化为平均量表达的形式，方程组方可封闭。常用的标准七一模型具有如下的 形式： 一万。z ，。婴一妻矿 一”一，i 可+ i j 一一v 2 - 2 v , 警一 矿 一b u 。旦盟 仃口8 工 一面；l 盟 8 其中t 。c k2 8 是湍流粘性系数，是湍流普朗特数。 将上面各式代入控制微分方程组并化简，可将每个方程均化为二维对流扩散方 程的标准形式： 判+ 趔。旦f r 盟1 + 三f r 丝1 + s ( 2 - 4 ) o x a y o xi o xj 砂l 妙j 各方程中，毋，r ，s 的表达式列于表2 - 1 中。 驴 r u 卢+ 卢。 一面o p + 磊o 【t 、p t 一、，i o uj + 专( 。+ 肛，) 警) 一气 y 芦+ p ， 一号+ 去( 。+ p i ) 詈) + 专( ( 一+ “) 孑) 一气 o p 。 0 _ 。 盯c 七 “+ 丝 g t p 1 - p p “ ，o r k s “+ 盟 e ，、c 。瓯一c 2 , 胪) + 僻 口f托 表2 - 1f e t l o c k s 模型各方程系数和源项 表中，p ，。p c ，譬，g t 。 2 ( 尝) 2 + ( 芳) 2 】+ ( 芳+ 尝) 2 。 对流扩散方程( 2 4 ) 的数值求解方法及控制方程组( 2 3 ) 的求解过程将在 下一章中介绍。 2 2 计算区域及边界条件 9 本文采用的计算域如图2 - 1 所示。 _z， b 制 图2 - 1 污水岸边排放问题计算区域 图中，河宽a b 为l ，计算流域长b c 为4 l ，e f 为排放口，宽度6 一面3 工，排 放口到上游起始边距离，- 嘉工。 边界条件如下： a b 边为入流边界条件，其上：u - u 。，( 入流速度) v - 0 ， 0 0 ， k 一盘。一0 0 0 3 7 5 u 0 2 ， 。箍。 b c 边上，e f 为喷口边界条件， u 一0 ，v 一圪，( 排放速度) 一0 。，( 排放口污染物组分浓度) k = k h 一0 0 0 3 7 5 e 。， 。慧。 1 0 其他部分及a d 边为固壁边界条件，其上对浓度使用导数边界条件墨o ， d v 法向速度矿取0 ，其他变量的固壁边界条件在湍流问题中需用壁函数法或低雷 诺数紊流模型处理。本文采用壁函数法，将边界条件直接引入到近边界节点的 离散格式中，详见第三章。 c d 边为出流边界条件，其上所有变量对流向坐标x 的导数为o 。 3 方程组求解的数值方法 上一章在紊流传质基本方程组的基础上引入了水深方向的平均化处理和河 底的阻滞因素，建立了污水岸边排放问题的控制方程组，并将各个方程整理为 形如( 2 4 ) 的二维对流扩散方程的标准形式。本章将讨论控制方程组的数值处 理方法。 对于粘性流体紊流模型的数值求解可以使用有限差分法，比如对k ，方程 采用r o e 格式进行离散，对其他方程采用一般方法离散，两个模块耦合迭代求 解。而对于本文讨论的此类问题，目前更为广泛采用的是有限体积法。与有限 差分法一样，有限体积法也是将微分方程通过一些手段离散到各个节点，形成 描述变量在相临几个节点上取值关系的代数方程，然后通过求解代数方程组得 到各离散点上待求变量的值。 要对实际紊流问题复杂的控制微分方程组应用此类方法，首先要作的处理 是解耦和线性化，即将每个未知量的方程写成形如( 2 4 ) 的标准形式，其系数 中含有的其他未知量及此未知量本身都作为常数看待。作这样的处理原因在于 目前的离散方法无法直接处理实际流体力学方程这样复杂的情况，但是对于那 些经典的模型方程( 比如这里的标准二维对流扩散方程) 的离散却是在理论上 有较多研究和保障的。而通过迭代，当方程中显式出现的变量与其在系数中及 其他方程中的取值很接近的时候，原复杂微分方程也就近似得到了满足。解耦 和线性化及相应的迭代过程在实际粘性流体数值模拟的一大类问题中被广泛用 到，在3 1 中将对此处理过程进行详述。 各方程经过解耦和线化处理，统一为方程( 2 4 ) 的形式后，采用何种方法 离散此方程是影响整个数值模拟精度的重要问题。显然，对此方程进行精度高 而数值性的非物理效应小的离散是研究不同湍流模型对流场细微作用的差异的 前提。有限差分法，有限体积法，以及近年发展起来的有限分析法等都可以对 此问题进行离散。但这些方法中，有限体积法得出的离散方程，要求因变量的 积分守恒对任意一组控制体积都得到满足，因此不论网格粗细，对攘个求解区 域，积分守恒也能得到满足，相比起来，有限差分等其他数值方法只能在网格 很细密时，积分守恒才能满足。综合比较起来，有限体积法有其优势，故为目 前同类问题使用得最为广泛的离散化方法。本文采用有限体积法对方程( 2 4 ) 进行离散，具体在3 2 中进行介绍。 另外一个问题是压力p 的处理，可看到在控制方程组( 2 3 ) 中，u ，v ， o ，k ，p 均有各自的方程，而密度p 和压力p 却没有。实际上，在考虑了能 量守恒方程的情形下，密度通过状态方程与温度相联系，在考虑组分浓度方程 的情形下，密度通过布辛涅斯克假设与组分浓度相联系。并且按照布辛涅斯克 假设，除动量方程中的体力项之外，密度可近似视为常数。压力p 与速度u ， v 联立的方程包括动量方程和连续性方程三个方程，但连续性方程中不显含压 力，故不能对压力直接实施3 1 中的迭代过程。也无法像其他一些问题中那样 采用某些特殊的方法使p 与u ，v 解耦，本文的问题中对求解压力需使用特殊 的方法，即p a t a n k a r 和s p a l d i n g 提出的s i m p l e 算法。相关内容在3 _ 3 中简介。 3 1 数值计算中方程组的解耦和线性化 如第二章所述，控制微分方程组( 2 。3 ) 中五个方程均可化为二维对流扩散 方程标准形式( 2 4 ) 。各方程中r ，s 的表达式如表2 - 1 所示，其中含有本方程 及其他方程的因变量。由于这种耦合和非线性，使方程难以直接进行离散。解 耦和线性化是指求解各个方程时将系数r ，源项s 视为相对不变的量，事先由， 备因变量初始假定分布或迭代当前分布直接求出。这样每个方程便均可按( 2 4 ) 的形式进行离散求解了。当然这种分步降低方程复杂性的方式必然需要迭代。 具体的迭代计算过程如下，首先对各变量假定一个初始分布，然后按u ，v ，o ， k ，s 的顺序计算各变量的方程的系数r 和源项s 求解各方程，并将解得的值 进行松弛处理后更新原先数据。反复求解各方程直到收敛。相应的流程图如图 3 1 。 图3 - 1 求解控制微分方程组的流程图 1 4 3 2 有限体积法求解二维对流扩散方程的细节 3 2 1 有限体积法简介 有限体积法( f i n i t ev o l u m em e t h o d ) 又称为控制体积法。其基本思路是： 将计算区域划分为一系列不重复的控制体积，并使每个网格点周围有一个控制 体积。将待解的微分方程对每一个控制体积积分，假定因变量值在网格点之间 的变化规律之后，对体积分进行数值积分，可得出一组离散方程。其中的未知 数是网格点上的因变量的数值，某些项通过g a u s s 散度定理，体积分化为控制 体积表面的面积分，这样也在离散方程中产生因变量在控制体积各边界节点上 的取值。经过上述过程，原来的偏微分方程组转变为各个网格点上的变量取值 间的代数方程组。求解代数方程组，便可求得在各个控制体积的体积分意义上 近似满足原偏微分方程组的问题的离散解。 有限体积法的基本思路易于理解，并且有限体积法的离散方程保持有明显 的物理意义：因变量在有限大小的控制体积中的守恒。这如同微分方程表示因 变量在无限小的控制体积中的守恒原理一样。有限体积法得出的离散方程，要 求因变量的积分守恒对任意一组控制体积都得到满足，对整个计算区域，自然 也得到满足。这是有限体积法最大的优点。有一些离散方法，例如有限差分法， 仅当网格极其细密时，离散方程才满足积分守恒；而有限体积法即使在粗网格 情况下，也显示出准确的积分守恒。 从积分区域的选取方法看来，有限体积法类似加权剩余法中的子区域法。 从未知解的近似方法看来，有限体积法属于采用局部近似的离散方法。就离散 方法而言，有限体积法可视作有限单元法和有限差分法的中间物。有限单元法 必须假定值在网格点之间的变化规律( 即插值函数) ，并将其作为近似解。有限 差分法只考虑网格点上的数值而不考虑值在网格点之间如何变化。有限体积法 只寻求结点值，这与有限差分法相类似；但有限体积法在寻求控制体积的积分 时，必须假定值在网格点之间的分布，这又与有限单元法相类似。在有限体积 法中，虽然使用到插值函数，但插值函数只用于计算控制体积的积分，得出离 救方程之后，便不再使用。并且如果需要的话，可以对微分方程中不同的项采 取不同的插值函数。 3 2 2 网格系统 有限体积法有自身的离散网格系统。其构造方式是把整个计算区域划分为 多个互不重叠的区间( g r i d ) ，然后确定每个区间上的节点位置以及该节点所代 表的控制体积。离散化之后，存在以下几种几何要素： 节点：最终求得的未知量取值的几何位置； 控制体积；有限体积法中将对原微分方程积分的范围： 界面：各个控制体积的分界面； 网格线：连接各个节点的( 曲) 线簇。 对于二维问题，标准的有限体积法结构化网格如图3 2 。 ，、 - ：、：、 t 矿 一p 童f 苴 9 ，强i f 二 j z石 图3 2 有限体积法离散网格系统 图中p ，s ，为节点，阴影部分为节点尸的控制体积， ，s ，e ， w 处为控制体积的界面。控制体积在两个方向上的宽度分别为缸，a y 。 在我们的问题中，虽然方程组( 2 - 3 ) 中每个方程都按图3 2 的网格系统进 行离散，且网格尺度一样，但是速度与其他标量方程的网格却不重合，而是在x 。 y 两个方向上均错开半个控制体尺度。其原因是如果速度与压强的网格重合， 则在p 点计算动量方程的压差项时，只能反映e ，处压强的情况，而完全使 用不到p 点处的值，从而压强在p 点处可能偏离e ，舻点的情况完全无法影响 速度。使用上述速度与其他标量交错的网格( 如图3 3 所示) 则可避免此情况。 臻 醚 瓣 醚露 鬻 誉 辫 ， 图3 3 各变量控制体积的相对位置 3 2 3 离散格式的构造 本节将介绍按有限体积法导出方程( 2 。4 ) 在图3 2 所示的网格下离散格式 的过程。 将方程( 2 4 ) 在图3 2 所示的控制体积内积分，有： 止掣+ 掣l d v c - 丘降( r 詈) + 专( r 耕e + 胖t ， 式中，积分体积记为p ，以便与速度y 相区分。 对方程左边的对流项使用g a u s s 散度定理( fv a d v 肠) ，化为控o n d s ja”舢s 制体界面上的面积分，并由物理量在界面中心上的值近似计算，有： ( p u c k ) , 一( j d u 妒) 。+ ( 卢l ，别l 一( 户l ，鲫) ， 式中，a e 。a 。，4 。，a ，为控制体右，左，上，下界面的面积。而角标e ， w ，撑，s 对于p ，u ，v ，矿表示其在相应界面上的值，由于最终目的是要形 成各个离散点间的方程，故界面处的值需通过插值化为节点值进行计算。有限 体积法中对流项的插值可用多种方式，考虑到我们求解的问题存在着强烈的主 流，故采用迎风格式有利于突出上游的影响，此时有： 丸- 知，屯一办。( u ，0 ) 妒，一妒p ，啦一九。( u 0 ) 九一妒p ，丸一九。( v ，0 ) 九一妒，丸一如a ( v t 0 ) 方程右边第一项( 扩散项) 作同样处理，化为： ( r 警一) 。一( r 署一) 。+ ( r 詈月) 。一( r 詈一) ， 扩散项的中的导数在有限体积法中贯采用中心差分离散到节点。离散后 有： m 锗州。锗姒饼啾锗 方程右边第二项( 源项) 本来在( 2 4 ) 的形式上是常数项，其积分就为鼬k ， 但一般为了利于格式的收敛性，都要作所谓源项线性化处理，即s = s 。+ s ，弗， 积分后为s c k + s p 蚱v o 。 经过以上处理，方程( 3 1 ) 化为节点间的代数方程： 口，九一口矿_ + 4 九+ 口5 九+ 拉v 妒+ 6 ( 3 2 ) 其中： ”簖批u ) 以) ”簖州p u 川 器x 帆) ”簖x m 帆) 口，- 口+ 口：+ 口，+ 口。+ ( p u ) 。4 一( p c ，) 。a ，+ ( p y ) _ 。一( p 矿) ，a ，一s ，t b - s c 屹 ( 3 2 ) 式即为标准形式的二维对流扩散方程( 2 4 ) 采用有限体积法得到的 离散格式。( 其中，对流项的插值采用迎风格式，且包含了u ，0 ，uc0 ，v 0 ， v t0 各种情况。) 1 8 3 2 4 边界条件的引入 将离散格式( 3 2 ) 对求解域中各节点写出，求解形成的方程组便可得到因 变量在全域节点上的值，但需要先对邻近边界的节点引入边界条件，对固定壁 面的情形，有限体积法中引入壁函数处理的方法是比较特殊的，需特别注意。 如2 2 所述，本文所求问题的边界条件包含入流及喷口边界条件，出流边界 条件，固壁边界条件三类。下面分别介绍其引入到边界点离散格式中的方法。 入流及喷口边界条件中，待求量在边界上的取值给定，将其直接引入边界 点离散格式即可，不需作其他处理。比如在流域左边界上，毋值已知，则对于 邻近左边界的节点，格式( 3 2 ) 中的如用边界值代替即可。除了入流边和排 放口，上下固壁边界上法向速度y o ，也按此方式直接引入。 出流边界条件要求备变量在边界上的法向导数为o ，这对于一次性直接求 解全域的情况，可能需要考虑虚拟节点法等其他方式，但我们求解全域的算法 是逐行计算再迭代( 见3 2 4 ) ，故对于法向导数为0 的条件的处理，可在每一 个迭代步中，将出流边界上的当前值像入流及喷口条件一样直接引入紧靠边界 的一列节点离散格式中解得此列节点取值后，再将边界值赋为与其相同，这 样当变量在紧靠边界的一列节点上收敛时，流向导数为0 的边界条件便得到满 3 凸 足。浓度在上下固壁处的边界条件! 兰一0 也按此法处理。 砂 【，k ，s 在上下固壁处的边界条件需使用壁函数法处理，且引入有限体积 法离散格式的方式是直接修正源项( 边界节点的联系被解除) ，表3 1 列出了固 壁影响对各个变量离散格式的源项修正量。 变量 a s ， a s c 一牮1 1 41 2 屯 u u “ o p c ，3 1 4 矿 u + 矿 监旷 七 y p 缈， ， 譬34 3 2 匕 0 由， 。 表3 - 1 固璧影响对各变量的源项修正 1 9 表中，u + 。三l n b + ) ，膏为i i ( , e r i t i e n 常数，对于光滑壁面有r 0 4 ，e 是 与表面粗糙度有关的常数，取9 8 。y + 为一无量纲距离，定义为生，f 尘，缈 y vp 为到壁面的距离，这里y + 是指节点p 处的y + ，网格尺度应取合适的值，使y + 大 于1 1 6 3 ，方可保证p 位于对数律层。p 点处y + 的计算式文献【1 1 】推荐采用公式 y + 鱼尘曼型，y ，为p 到壁面的距离。 岸 爿。为控制体积在壁面处的面积，在本文的问题中为缸。七；为尸点前次迭 代的| 值。 一p q “4 k1 2 u ，d + 为壁面剪切应力。 这样，u ，七，近边壁节点的离散格式引入固壁边界条件后变为如下形式： ( 口，一磷k 砩一4 ，如+ 九+ 如+ 6 + 峨( 上边壁) ( 口，一船p 屹k 一口如+ 九+ 口“+ 6 + 丛c ( 下边壁) 3 2 5 全域节点变量的求解 格式( 3 - 2 ) 为五点格式，求解时需将整个求解域中所有内节点的离散方程 联立方可封闭。但这样形成的代数方程组是五对角的，求解较麻烦。按下述逐 行迭代的方法，可将运算化为只涉及求解三对角方程组，实施起来更简单。逐 行迭代是指对全场所有内节点上的妒设定初始值，然后将格式( 3 2 ) 中如，九 视为已知，用妒的初始分布或前次迭代值代替，这样格式简化为三点格式，在 一列节点上就可以求解。现从左到右，逐列求解各节点的值，但注意到此时节 点间横向的联系，只是从左到右得以实施( 求解一列节点时，如是求解过的， 传递着最左边边界条件的影响，而以是初始随意设定的) ，故需要反复迭代求 解，即从左到右逐列求解的过程需反复实旌。在迭代的过程中，右边界条件的 影响逐渐传向左方。当变量在各列的取值均收敛后，便完成了整个区域的求解。 按上述方法，求解每列节点时，都需解形如下式的三对角代数方程组： 办- c 1 一声2 疵+ 见丸- o f 2 丸一c 2 一只办1 + d i 办一口旃“- c i 一成丸。1 + d j 丸一口。丸+ l - c 以+ 。- q + 。 ( 此形式对应3 2 3 中介绍的各种边界条件，因为网壁边界条件虽然解除了 哦和九+ 。的联系，增加了源项的修正，但也可用形式上的蛾和九+ 。来实施源项的 增量。) 上式中除了第一个和最后一个方程外，其余方程都具有以下形式： 一成珐- 1 + d i 唬一a 商+ 。一c i