（热能工程专业论文）二维开口方腔混合对流换热的数值模拟.pdf

哈尔滨工程大学硕士学位论文 摘要 开口方腔的对流换热问题是一基本的物理现象，工程技术上的许多问题 都可以简化成这一基本问题的求解。本文从探讨这一问题的物理原理出发， 推导其基本控制方程，为了更深刻的分析计算结果本文引入了相似原理，针 对几种基本物理模型，通过改变雷诺数、格拉晓夫数及纵横比系数进行了大 量的数值计算，通过对流线、等温线及方腔开口面温度、速度曲线的分析， 深刻的理解了这一物理现象，最终得到了控制流场及温度场的最主要参数是 理查德数。 关键词： 开口方腔；混合对流：数值模拟 哈尔滨t 程大学硕士学位论文 a b s t r a c t m i x e d c o n v e c t i o nh e a tt r a n s f e ri no p e n e d e n d e dc a v i t yi sab a s i c a l l y p h y s i c a lp h e n o m e n o n m a n yt e c h n i c a lp r o b l e m sc a nb es o l v e dt h r o u g hs t u d y i n g t h i sb a s i cp r o b l e m i nt h i s p a p e r ，t h ep h y s i c a lp r i n c i p a li s s t u d i e da n dt h e d o m i n a t e de q u a t i o ma r ed e d u c e d i no r d e rt os t u d yt h i sp r o b l e mp r o f o u n d l y ，t h e s i m i l a rp r i n c i p a li s i n t r o d u c e d b yc h a n g i n gt h er e y n o l d sn u m b e r ，g r a s h o f n u m b e ra n dt h ea s p e c tn u m b e r , ag r e a td e a ls i m u l a t i o na l em a d ef o rs o m e b a s i c a l l yp h y s i c a lm o d e l s b ys t u d y i n gt h es t r e a m l i n e s ，t h ei s o t h e r m sa n dt h e t e m p e r a t u r ea n dv e l o c i t yp r o f i l e sa tt h ea p e r t u r ep l a n eo ft h eo p e n - e n d e dc a v i t y t h ec o n c l u s i o ni sg a i n e dt h a tt h ep a r a m e t e rt h a tc o n t r o l st h ei s o t h e r m sa n dt h e s t r e a m l i n e sm o s t l yi sr i c h a r dn u m b e r k e y w o r d o p e n - e n d e dc a v i t y ；m i r e dc o n v e c t i o n ；n u m e r i c a ls i m u l a t i o n 哈尔滨工程大学 学位论文原创性声明 本人郑重声明：本论文的所有工作，是在导师的指 导下，由作者本人独立完成的。有关观点、方法、数据 和文献的引用已在文中指出，并与参考文献相对应。除 文中己注明引用的内容外，本论文不包括任何其他个人 或集体已经公开发表的作品成果。对本文的研究做出重 要贡献的个人和集体，均已在文中以明确方式标明。本 人完全意识到本声明的法律结果由本人承担。 作者( 签字) ：b 垡趁 日期：日1 年 月勺日 哈尔滨工程大学硕士学位论文 第1 章绪论 1 1 研究对流换热的必要性与重要性 在现代科学技术领域里，无论动力、冶金、石油、化工、材料、制冷 等工程，还是空间、电子、核能等高技术，都不可避免的涉及到加热、冷 却和热量传递的问题。在能源的开发利用中，热量传递现象更为普遍。如 何改善、提高热传播的速率，以达到使用最经济的设备来传递规定的热量， 或是用最高的热效率来实现能源合理利用的目的。 由于生产和科学技术发展的需要，这一方面的研究工作获得了广泛的 重视和发展。 首先，世界面临着能源短缺的局面。开发新能源和节约能量消耗引起 了世界各国有关部门的普遍关注。设计和制造各类高性能换热设备是经济 的开发和利用能源的最重要手段，这对于动力、石油、冶金、化工、制冷 及食品等工业部门都有着极为重要的意义。 其次，随着航空、航天及核聚变等尖端技术的发展，各种热设备的工 作温度也在不断提高，为了保证热设备能有足够长的工作寿命和在高效率 下安全运行，必须可靠而经济的解决高温设备的冷却问题。 最后，随着巨型计算机的迅速发展，密集布置的大功率电子元件在电 子设备中的释能密度日益提高。电子元件的有效冷却，是电子设备性能和 工作寿命的必要保证。因此强化传热已成为保证如大功率微电子元件等强 热流部件的可靠性以及实现如低品位热能开发利用等小温差有效传热的技 术关键。 1 2 国内外的研究概况 混合对流换热问题可以归结为浮力与惯性力共同作用下的对流换热问 题。对于这一问题近几十年来国内外已经进行了很多的研究”，但大都 哈尔滨工程大学硕士学位论文 集中在封闭方腔内的对流换热问题，研究浮升力与重力共同作用下的混合 对流换热问题。而对浮升力与其它惯性力共同作用下的混合对流换热问题 研究的较少，国内很少见到这方面的文献，而国外近几十年来针对开1 3 方 腔的混合对流换热问题做了不少研究。p e n o t 推导了二维开口方腔内混合 对流换热的控制方程并进行了数值求解，他的研究结果表明当方腔内部的 格拉晓夫数超过1 0 5 时就会产生不稳定的解；c h a n 和t i e n 等人集中研究了 瑞利数在1 0 3 1 0 9 之间的层流自由流动换热问题：关于开口方腔自然对流 的实验研究主要以h u m p h r e y 、s e r n a s 、k y r i a k i d e s 、s h o w o l e 以及t a r a s u k 为代表2 2 m 1 。关于这一问题更广泛的研究是由k h a l i lk h a n a f e r 、k a m b i z v a f a i 、e t t e f a g h 以及s o n i ab e s b e s 进行的。他们研究了开口方腔内部流场 特性，并研究了扩展计算区域对开口方腔内部的流场及换热效果的影响， 他们的研究结果表明数值计算时外延计算区域需要很大，同时他们的结果 也表明远处流场流体特性与使用的边界条件关系密切 1 2 - 1 5 】。m h i r i 等人针 对具有竖直射流的开口方腔混合对流换热问题进行了数值模拟研究，计算 结果表明方腔中流体的进出可以通过加大射流空气速度很好的避免，从而 达到很好的隔热保温效果。后来m h i r i 等人使用有限差分方法又对水平射 流的情况进行了研究，计算结果同样表明可以通过增大射流速度达到很好 的保温效果。 1 3 研究内容 该论文主要通过使用商用软件p h o e n i c s 对开口方腔的混合对流换热问 题进行研究，通过与参考文献【1 2 - 15 1 中的经典算例结果进行比较验证使用软 件的可行性，然后主要针对下面三种情况进行研究：( 1 ) 水平开口方腔的 混合对流换热问题，其中包括水平射流和外界气流水平流经方腔开1 3 的情 形。( 2 ) 单层竖直开口方腔外界气流不同方向角度流经方腔开口的混合对 流换热问题。( 3 ) 多层开1 3 方腔的对流换热问题，研究不同层之i l 丑j 的相互 影响。通过改变相似参数，对不同的情况进行比较最终提出合理的加强对 流换热和保温隔热的方案。 2 第2 章物理现象的描述及控制方程 2 1 基本物理概念介绍 2 1 1 对流换热 对流是指由于流体的宏观运动，从而流体各部分之间发生相对位移、 冷热流体相互掺混所引起的热量传递过程。对流只能发生在流体中，而且 由于流体中的分子同时进行着不规则的热运动，因而对流必然伴随有导热 现象，我们称流体流过一个物体表面时的热量传递过程为对流换热。 就引起流动的原因而论，对流换热可区分为自然对流与强制对流两大 类。不依靠泵或风机等外力推动，由流体自身温度场的不均匀所引起的流 动称为自然对流。例如，暖气管道的散热、不用风扇强制冷却的电器元件 的散热以及事故条件下核反应堆的散热都是自然对流换热的应用实例。不 均匀温度场造成了不均匀密度场，由此产生的浮升力成为运动的动力。如 果流体的流动是由于水泵、风机或其他压差作用所造成的，则称为强制对 流。冷油器、冷凝器等管内冷却水的流动都由水泵驱动，它们都属于强制 对流。 在强制对流换热中有时需要既考虑强制对流亦考虑自然对流，应用相 似分析法可知，g ，数中包含着负升力与粘滞力的比值，而由惯性力与粘滞 力对比可得如数，我们需要的是浮升力与惯性力的对比，引入理查德数胁 ( r i = g r r e 2 ) 。一般认为，r 扛g ，r e 2 0 1 时自然对流的影响不能忽略， 而尉= g r r e 2 1 0 时强制对流的影响相对于自然对流可以忽略不计【1 6 】。 2 1 2 浮力射流 一股流体从几何尺寸远小于受纳流体所占空间尺寸的喷口流入受纳 流体，并同其混合得流动状态，叫做射流。射流有多种不同的分类方法， 哈尔滨工程大学硕士学位论文 而我们文中主要关心的是无限静止环境中不可压缩湍射流、湍卷流和湍浮 力射流。纯射流也称动量射流，从喷口流出的流体与其周围的受纳流体物 理性质相同，射流的形成主要是喷口处初始动量的作用。卷流也称羽流， 射流的形成是喷口处流体与受纳流体的密度差产生浮力的作用。浮力射流 也称强迫卷流，射流不仅有初始动量，而且还受浮力的作用1 1 7 】。 | 厂、 、” |鞭m ，一j 警拄! 飞1 a ，一 五 c a ) 纯射流( b ) 卷流 ( c ) 浮力射流( d ) 负浮力射流 图1 - 1 静止均匀环境中的射流 浮力射流的运动既受初始动量的影响，又受浮力的作用。显然纯射流 和卷流的是它的两种极端情况。本文主要关心与讨论的是均匀环境中的浮 力射流。均匀环境指的是，环境流体的密度各处都是相同的。设有一股密 度为风、初速度为的流体，从喷口铅直 的射入密度为见的静止均匀环境中。在喷 口处单位体积流体的动量为p o u o ，所受的 浮力为( 见一p o ) g ，g 为重力加速度。若 见 p o ，则浮力的方向铅直向上，叫做正 浮力。在正浮力作用下的射流，就叫做正 浮力射流。若成 0 ，即射流温度大于环境温度，方程中最后一项是正浮力作用； 反之为负浮力作用。 ( 3 ) 能量方程： 在满足不可压缩性条件下，略去粘性耗散项，能量方程( 2 - 3 ) 可写为 塑+ u 塑：a 旦i 旦l + 旦 ( 2 一7 ) + 一2 a l l + l z j a i 瓠ie ) j 瓠i ) cp 式中a ：l ，叫做流体热扩散系数，这就是布辛涅斯克型能量方程， p 口o p 现在我们来讨论能量方程粘性耗散项可以略去的条件。此项的量级是 2 艺一p 等，而在运动方程中，浮力项与粘性力项是同量级的，其量级是 百u p , , g f l a t ，因而粘性耗散项与能量对流项之比的量级是： 墨丝丝墨丝 p 乒p u a t lc p 若墨丝1 ，则在能量方程中可略去粘性耗散项，对于空气中的浮力射流， 这一条件一般是满足的。 通过上面讨论得到的( 2 4 ) 、( 2 6 ) 及( 2 7 ) 就是引入布辛涅斯克假 设后得到的控制方程，现在单列出如下： o u l ：o o x 。 8 ( 2 4 ) 凡( 警m 考 _ 一瓦o p + 蕞+ 舰酗丁 c z 勘 o 西t + u 飞0 7 = 口考+ 号 协， 因此二维、不可压缩流体的常物性无源稳态流动 换热的方程为： 质量方程：塑+ 堡：0 融砂 直角坐标系下对流 动量方程：z ，瓦o u + v 万o u = 吉p 挈o x + v 窘+ 乎o v 册 删l 血+j & o x v 堡o y p o x 豺o y 酬)印。l 2。2j ”1 “， 能量楹“一s t + v a t o x o y 窑o x 芬 【 2 却2j 2 3 控制方程的无量纲化 引入下面无量纲量 x ：三 日 臼：堡量， t 一2 】，：上， n p 砰 侧 u ：士， ( a h ) v ： ! h ) ( 2 8 ) ( 2 9 ) ( 2 1 0 ) ( 2 1 1 ) 可以得到无量纲控制方程： 质量方程：o e _ u + 关：o ( 2 - 1 2 ) o xa y 动量方程：u 型o x + 矿丝一面o p + p r i ( 描0 2 u 。叶焉0 2 u o yo y ) ( 2 - 1 3 ) 捌l 描2 2j u 筹+ 矿万o v 一万o p + 脚( 豢+ 祟o y 卜r 衍p c z 出， 硝a ya yi 硝2 2j 。 一。4 9 哈尔滨工程大学硕士学位论文 艟方程：u 筹+ 矿石0 0 = l 硝0 2 0 i + 萨a 2 0 协 式中的无量纲相似参数 r a = g r p r , 西：型三：堡坦 v 2 4 本章小结 j p ，：一v ， r p o r 本章包含以下几个方面的内容： ( 1 ) 、自然对流换热、强制对流换热及混合对流换热的基本概念，并 引入了理查德数，对自然对流换热与强制对流换热强弱判定方法进行了讨 论； ( 2 ) 、射流的定义以及分类，并着重介绍了浮力射流； ( 3 ) 、相似原理的简单介绍及相似参数的引入。 ( 4 ) 、浮力射流的基本控制方程，并推导了引入布辛涅斯克假设后的 控制方程，最终给出了直角坐标系下二维、不可压缩流体的常物性无源稳 态流动的对流换热的控制方程。 ( 5 ) 、引入无量纲数，对控制方程进行无量纲化。 1 0 哈尔滨工程大学硕士学位论文 第3 章开口方腔对流换热的数值解法 计算流体力学和传热学自七十年代以来的成就，显示出它在人类深入 研究流动与传热现象，以及在工程应用方面的强大生命力。在航空、气象、 海洋、流体机械、建筑和车辆设计等各个领域显示出巨大威力。近年来， 由于计算机速度和储存信息能力的大幅度提高，由于计算机自动生成三维 物体网格能力的迅速发展，计算机软件水平突飞猛进，给科学发展和工程 应用设计带来了根本性的变化。 计算流体力学与传热学作为一门崭新的学科，它主要用来研究和描述 各种流动与传热现象，其目标是在工程上用数值试验代替实物实验，用计 算机模拟和设计自然环境、生物体和工程机械装置。计算机的大量使用， 计算流体力学和传热学通用大型软件的成熟，彻底改变了人们对于实验和 设计的传统观念。通过数值模拟对工作过程细节的了解，使得工程装置的 优化设计现在已经作为一种新的手段，可以预报真实的流体机械、换热与 燃烧装置、工业炉、大气污染等现象的全过程。可以得到设计所需的各种 定量数据，又能把实验所需的人力、物力及财力减到最低限度，实现了真 币意义的设计革命。通用计算机软件的发展，使数值模拟成为与模拟实验 同样可靠和有效的手段。在计算机程序中，改变几何、物理和化学因素， 以及改变环境条件都比在实验室中改变这些参数要容易得多。既可以加快 研究速度，又可以拓宽研究参数的变化范围，从而增加科研的深度与广度。 p h o e n i c s ( p a r a b o l i c ，h y p e r b o li c o y e l l i p t i c n u m e r i c a l i n t e g r a t i o i lc o d es y s t e m ) 程序是英国c h a m ( c o n c e n t r a t i o nh e a ta n d m o m e n t u ml t d ) 公司完成的享誉世界的通用计算流体力学与传热学软件， 是英国帝国理工学院著名学者d b s p a l d i n g 教授及4 0 多位博士集2 0 多年 心血的典范工作。它是流体力学、气体动力学、传热传质、化学反应及相 关学科的通用计算软件。可以计算二维、三维、稳态、瞬态、粘性、非粘 性、可压、不可压、层流、湍流、单相及多相流等各种流动与传热现象。 目前，该程序已广泛应用于航空航天、船舶、水利、汽车、环保、发电、 哈尔滨工程人学硕士学位论文 冶金、核工业、生物医学和农业等行业。 本文就是运用p h o e n i c s 通用程序为平台，通过编写程序求得丌口方 腔混合对流换热的数值解。 3 1 控制方程 3 1 1 控制方程的推导 对于一个三维的对流换热过程，控制方程有三种即质量守恒方程、动 量守恒方程、能量守恒方程。质量守恒方程和能量守恒方程很容易根据质 量守恒定律和能量守恒定律获得，在这里我们直接引出如下： 一、质量守恒方程： 竺+ 掣+ 掣+ 掣：o ( 3 - 1 ) o t缸加瑟 也可用矢量符号写出为： 詈+ a i v ( p v ) - o ( 3 - 2 ) 对于不可压缩流体，其流体密度为常数，质量守恒方程简化为： 咖( v ) = 0 ( 3 - 3 ) 二、温度形式的能量守恒方程： 掣砌( 妒脚v 肿卜 。， 对不可压缩流体有： 鲁+ 柳( w ) = 柳【去矿们r j + 軎 一s ， 三、动量守恒方程： 动量守恒方程也就是纳维尔一斯托克斯方程，为了便于推导这里先介绍 一些基本概念： ( 1 ) 、速度增量： 1 2 哈尔滨工程大学硕士学位论文 阮= 瓦o u , 8 x ，：1l ( 虿o u i + 善卜- + 虱1 ( 虿o u , 一- - 警 1 6 x j = 8 x ：+ c o , s 啦 式中对称张量白= 虱1 ( 钙a u , - + 等 ；反对称张量= 丢( 考一鼍 ( 2 ) 、在粘性流体的一般运动中，通常可将应力张量分成两个部分： p o _ 一p 磊 + 厩 ( 3 - 6 ) 其中- p 4 ) 是代表流体静止时应力状态部分( p 是静压力，与热力学压力 相同) ； p ： 是由于粘性流体运动产生的应力状态部分，与变形速度有关 称它为偏应力张量。当流体静止时，偏应力张量为零。 ( 3 ) 、斯托克斯假设： ( a ) 粘性流体中应力是变形速度的线性齐次函数； ( b ) 流体是各向同性： ( c ) 流体静止时，流体中的应力就是流体静压力； ( d ) 平均法向应力与体变形速度无关。 由斯托克斯的前两个假设，可将偏应力张量表示如下 - 五 气 + 2 k ( 3 7 ) 式中五和是流体的特性参数，= e 1 1 + p 2 ：+ e 3 ，是体膨胀速度。由( 3 6 ) 式得总应力张量 p ， = ( 一p + 2 ) 磊 + 2 k ( 3 8 ) 上式是牛顿流体的本构方程( 表示应力与变形速度之间的关系式) ，也称为 广义牛顿摩擦定律。从( 3 7 ) 式求得偏应力的法向分量平均值 ；瓦= ( z + 弓 = 矗 c ，。， 哈尔滨_ 程大学硕士学位论文 其中系数凡( = a + 詈) 表示偏应力的法向分量平均值与体膨胀速度之比的 比例常数，叫做体积粘性系数( 或第二粘性系数) 。从( 3 - 8 ) 式求得法向 应力平均值 ；p i = 一p + ( 丑+ 詈 e 。( 3 - 1 0 ) 由上式可见，压力一般不是平均法向应力。按照斯托克斯的第四个假 设平均法向应力与体变形速度无关，则我们必须取 z + 弘。 因此：五= 一三“ 3 从而在可压缩流体中，系数 和动力粘性系数a 联系起来，所以应力 与变形速度之间关系式( 3 8 ) 变为： = 一詈触) ”z ( 3 - 1 1 ) 这个关系式只在体变形速度很小的情况下才是适合的，而对于大多数 流体来讲，体膨胀速度般不大，所以造成的误差不大。 对于不可压缩牛顿流体，此时由= 0 ，则本构方程( 3 - 1 1 ) 式可简 化为： 岛) = p 岛 + 2 ( 3 1 2 ) 根据牛顿第二定律，可写出直角坐标系的x ，y ，z 方向的运动方程式 掣= m + ( 誓+ 等+ 鲁 掣= 以十( 鲁+ 誓+ 誓 掣= 蹦+ ( 等+ 等+ 割 1 4 ( 3 】3 ) 詈+ 撕( w ) 2 六+ 咖( 憎a d “) 一吉塞7 考；+ 西v ( v v ) 2 + 西v ( 增r a d v ) 一吉考 c ，t 4 ， 詈瑚v ( w v ) = f + d i v ( 删小划 我们现在将质量守恒方程、动量守恒方程、能量守恒方程写在一起， 如下： 掣砌( 删训+ 咖( 脬砒) 警+ 瓯 旦a t - + 拆( p v v ) = 啊十旃v ( t g r a d v ) 宇+ e f 掣胁( 州= 疆砌( 炉d 砂一害+ 挈砌卟面v 陪d ，卜 ( 3 1 5 ) 其中= 孔窑一d i v v ) + 号( 割+ 孔割 e = 昙( 雾) + 号( p 考一d i v v l + 昙( 考) 凡= 丢( 巷) + 号( 害) + 丢( 西o w 一；西v v ，? 于粘性系数为常数的不可压缩流体，配；s v ：瓯：o ，( 3 1 5 ) 式可 简化为 。 哈尔滨工程大学硕士学位论文 d i v ( v ) = 0 罢+ d i v ( 胛) = 六十硪v ( y g r a d ”) 一土p 望c 2 r 詈+ 柳( 州) = + 旃v ( 卿d v ) 一吉考 警+ 机( w v ) 2 正胁( v g r a d w ) 一j l 瓦o p 詈+ 蹄v c v 丁，= 咖( 去g 豫d 丁 + 告 ( 3 1 6 ) 上面五个方程中包含了6 个未知量，“，v ，w ，p ，t 及p ，还需补充一个 联系p ，p 的状态方程，方程组才能封闭： p = f ( p ，t ) ( 3 - 1 7 ) 对理想气体可有： p = p r t ( 3 - 1 8 ) 其中r 为摩尔气体常数。 3 1 2 控制方程的通用形式 在流动与传热问题求解中所需求解主要变量( 速度及温度等) 的控制 方程都可以表示成以下通用形式： 堕d 盟+ v ( 矽) = r ( r g 船d ) + 只(319)t d id i 式中为通用变量，可以代表 ，v ，w , t 等求解变量；f 。为广义扩散系 数；咒为广义源项。这里引入“广义”二字，表示处在r 和墨位置上的项 不必是原来物理意义上的量，雨是数值计算模型方程中的一种定义，不同 求解变量之间的区别除了边界条件与初始条件外，就在与l 与咒的表达式 的不同。式( 3 1 9 ) 也包括了质量守恒方程，只要令妒= 1 ，i _ ；= 0 ，s 。= 0 即可，下表给出了可压流场控制方程的，f 。及& 值。 1 6 哈尔滨t 程大学硕士学位论文 表3 - 1 可压流场控制方程的矿，f d 及足值 对于粘性系数为常数的不可压缩流体= 鼠= s w = 0 控制方程有两种形式即守恒型控制方程和非守恒型控制方程。式 ( 3 - 1 9 ) 是通用控制方程的守恒型式，从数值计算的观点，守恒型方程有 两个优点。在计算可压缩流动时，守恒型的控制方程可以使激波的计算结 果光滑而且稳定，而应用非守恒型控制方程时激波的计算结果会在激波前 及后引起解的振荡，并导致错误的激波位置。所以空气动力学的数值计算 中守恒型的控制方程特别受到重视，并用通量的列矢量的一阶导数的方程 组的形式来表示。在计算传热学中，希望数值计算结果能满足守恒定律， 而要保证做到这一条，应该采用守恒型控制方程。也就是说只有守恒型的 控制方程才能保证对有限大小的控制容积内所研究的物理量的守恒定律仍 然得到满足。而且凡是从守恒型的控制方程出发，不论节点布置的疏密程 度如何，采用有限容积积分法导出的离散方程具有对任意大小容积守恒的 特性。离散方程的守恒特性是工程计算所需要的。 3 2 通用控制方程的离散 文献”1 中有详细的离散方程推导，我们这里直接引入二维问题的离散 化方程 d p 争p = n e 毋e + a w 纯+ q s 牵s + d n 牵n + b 3 - 2 0 式中： = 见爿( 吲) + i 只，0 1 = 仇一( 吲) e ，0 1 = 见4 ( 吲) + i 只，0 ( 缸) 。( 如) 。 = d 、彳( f 只f ) + 只，0 j图3 1 直角坐标系下的网格系统 口o 。：逊 。 a t b = 品缸缈+ 群 a p = d 8 + n w + a s + d n + 砖一s p 奴姆 e = ( 眺缈，d e = f e x ) 。 f = ( p u ) 。衄， e = ( p v ) 。缸 c = ( p v ) 。a x 贝克列数定义为 r ) 一f w a y 一 【万x ) 。 d a = r a x p n j b = r 5 y ) ， 1 8 妇) 妇) 哈尔滨工程大学硕士学位论文 只d ，只土d w ，只寺只= 每 【l 符号表示最大值。 函数一( h ) 可以按照所要采用的方案由下标选取 表3 - 2 各种不同方案( 格式) 的函数爿( 1 p 1 ) 互墨! 塑苎!翌兰i 蚓2 竺竺苎 中心差分 l 一0 5 1 j i 上风 混合 幂函数 指数( 精确解) 0 ，1 - 0 5 4 0 ，( ，刮p 州 i p l e x p ( 1 p 1 ) 一1 3 3 初始条件与边界条件 上面所讨论的控制方程及其离散形式适用于所有牛顿流体的流动与 换热过程，各个不同过程之间的区别是由初始条件及边界条件( 统称为单 值性条件) 来规定的。控制方程及相应的初始与边界条件的组合构成了对 一个物理过程的完整的数学描写。 初始条件是所研究现象在过程开始时刻的各个求解变量的空间分布， 必须予以给定。对于稳态问题不需要初始条件，所以本文对其不进行详细 的讨论。 边界条件是求解区域的边界上所求解的变量或其一阶导数随地点及 时间的变化规律。下面我们以图3 - 2 所示的突扩通道中有回流的流动为例， 来详细讨论边界条件的设定。一般的说，可以把计算区域的边界分成4 种 类型【1 9 】： ( 1 ) 、入口边界：入口边界上的条件必须给定，一般是规定了入口边 1 9 哈尔滨工程人学硕士学位论文 界上班的分布。 ( 2 ) 、中心线( 对称轴) ：在这 里有 v _ o 宴：o ，掣：o u u y ( 3 ) 、固体壁面：固体边界上又 可分一、二、三这三类边界条件。当 固体壁面为非渗透性时，壁面上“= v = 0 。图3 - 2 突扩通道 第三类边界条件上规定了边界上的庐值与娑( n 为法线) 之间的关系。当 册 计算区域中的流体与分隔壁外的流体有热交抉， 类型。 ( 4 ) 、出口边界：这是最难处理的边界 条件。目前广泛采用的一种处理方法是假定 出口截面上的节点对第一个内节点已无影 响，因此可以令边界节点对内节点的影响系 数为零。这样出口截面上的信息对内部节点 的计算就不起作用，也就不需要知道出口边 界上的值了。这种处理的物理实质相当于假 旦壁面很薄时就属于这一 出口边再 定出口截面上流动方向的坐标是局部单向的，图3 - 3 出口截面与内节点的关系 如图3 3 所示，与出口边界上e 点相邻接的第一个内节点p 与点之闻的关 系是通过尸点的系数盯。来规定的。如果对流作用比较强烈，则扩散作用可 以不计；又因为e 在尸的下游，e 对p 的影响亦可忽略。所以局部单向化 的假定导致口，= 0 。为了在数值计算中应用这一简化处理方法而又不致引 起过大误差，应做到： ( a ) 、在出口截面上无回流； ( b ) 、出口截面应离开感兴趣的计算区域比较远。 哈尔滨工程大学硕士学位论文 3 4 交叉网格 由前面对控制方程的推导、方程的离散及边界条件的介绍可以总结出 物理问题数值求解的基本过程，如下图所示： 线 性 闽 题 图3 - 4 物理问题数值求解的基本过程 由于工程中所遇到的流动与传热问题一般多为有回流情形椭圆型问 题。在求解这类问题时，可以用速度、压力( 或密度) 作为基本变量，也 可取涡量、流函数作为变量。前一类方法称为原始变量法，后一种方法叫 涡量函数法。本文使用的是原始变量法。 3 4 1 并置网格引起压力速度失耦 上面对控制方程的离散是在并置网格上进行的，并置网格是指把速度 “，v 及压力p ( 包括其他所有标量场及物性参数) 存储于同一网格上的网 格系统。纳维一斯托克斯方程是非线性的，它们的离散方程系数中包含有被 2 l 哈尔滨工程大学硕士学位论文 求变量“，v ，因而整个问题的数值求解必然带有迭代的性质。但从原则上 说，这并不构成特殊的困难。动量方程数值求解中所遇到的主要问题是与 高 一阶导数项兰的离散有关的： 啦 ( 1 ) 、采用常规的网格及中心差分来离散压力梯度项时，动量方程的 离散形式可能无法检测出不合理的压力场。 ( 2 ) 、压力的一阶导数以源项的形式出现在动量方程中。采用分离式 求解变量的离散方程，由于压力没有独立的方程，需要设计一种专门的方 法，以使在迭代求解过程中压力的值不断地得到改进。 以上两个关键问题都与压力梯度的离散及压力的求解有关，统称为压 力与速度的耦合问题。为了克服压力与速度间的失耦，可以采用交叉网格； 这样在采用分离式求解方法时各类变量能同步地加以改进，以提高收敛速 度。 3 4 2 交叉网格上速度分量位置的安排 所谓交叉网格就是指把速度“，v 及压力p ( 包括其他所有标量场及物 性参数) 分剐存储于三套不同网格上的网格系统。其中速度“存于压力控 制容积的东、西界面上，速度v 存于压力控制容积的南、北界面上，“，v 各 自的控制容积则是以速度所在的位置为中心的，如下图3 - 5 所示： ( a ) 主控制容积( b ) u 控制容积( c ) v 控制容积 图3 - 5 交义网格 哈尔浜上程大学硕士学位论文 由图3 5 可见，“控制容积与主控制容积( 即压力的控制容积) 之间 在x 方向有半个网格步长的错位，而v 控制容积与主控制容积之间在y 方向 有半个网格步长的错位。交叉网格这一名称即由此而来。 在交错网格系统中，关于“，v 的离散方程可通过对 ，v 各自的控制容积 作积分而得出。这时压力梯度的离散形式对为( 既一所) ( 出l ，对为 ( p n p e ) ( 却) 。，办即相邻两点间的压力差构成了譬，字，这就从根本上 o y 解决了采用一般网格系统时所遇到的困难，也是交叉网格成功的经验。 3 4 3 交叉网格上动量方程的离散 在交叉网格中，一般变量的离散过程及结果与并置网格中所述的一 样。但对动量方程而言，则带来一些新的特点，主要表现在以下两个方面： ( 1 ) 、积分用的控制容积不是主控制容积而是“，v 各自的控制容积。 ( 2 ) 、压力梯度项从源项中分离出来。例如对“。的控制容积。该项积 分为 i 夏一塞) 级方= 一! ( p 悟) 匆z ( p r - p 。) 匈 这里假设在“。的控制容积的东、西 界面上压力是各自均匀的，分别为p 。及 p ，。于是关于虬的离散方程便有以下形 式： 吼= a n b l 2 r i b + 6 + ( 西一p e ) a , 其中； b l m 是“。的邻点速度( 图3 - 6 中的”。，“。，u ，及“，) ；b 为不包括压力在 内的源项中的常数部分，对非稳态问题 为6 = & 矿+ 吒o o ；4 = 母是压力差的 作用面积；系数日。的计算公式取决于所 采取的格式，如前面章节所述。 ( 如) + 如) 一( 赫) + 卜 + 叫 f ( 卦) 1 4 - l ( 融仉i刎) 。十 ，一 - - i 彬1 1 缮 e t 妊l 图3 - 6 也的四个邻点 t 峄上 哈尔滨工程大学硕士学位论文 类似的，对心的控制容积作积分可得： = “。+ 6 + ( 肌一p ；) 4 3 4 4 交叉网格上的插值 在交叉网格上采用控制容积积分法来导出离散方程时，各控制容积界 面上的流量、物性参数等常需要通过插值的方式来确定。需要插值的量有 以下三类： ( i ) 、界面上的流量 例如控制容积的西界面上的流量耳可以按照心，甜。位置上的e ，r 插值而得： 张警+ r 譬如吐缈等+ ( p 0 。缈譬 而虬的北界面上的流量e 一。则可以看成分别由k 及在各自的流动 截面内的流量迭加而成： 只一。= ( ) 。( 撕) 。一+ ( ) 。( 8 x ) 。 以上两式中截面上的密度都要经过插值才能确定。 ( 2 ) 、界面上的密度可以采用线性插值法确定。如成可表示为： 见= 既储+ 所器 ( 3 ) 、界面上的扩散系数( 或扩导) 。如“。的北界面上的扩导或一。可 以表示成：址= 器+ 黯= 盘+ 盘 r 。 r h ! rpr n r e i v ： 其中r 。，r 。，r p 及r 。都是节点上的扩散系数。 3 4 5s i m p l e 算法 在获得动量方程的离散形式后，如果采用“，v ，p 同时求解的方法，则 2 4 哈尔滨j ：程大学硕士学位论文 需将连续性方程在主控制体上离散( 采用交错网格时速度均在主控制的界 面上，不必插值) ，然后用直接解法计算在给定的一组系数下各节点“，v ，p 之值。根据计算所得的新值，改进代数方程的系数，再用直接法解出与新 的系数相应的”，v ，p 值。如此反复，直到收敛。这种“，u p 在全域范围内同 时求解的方法由于要耗费巨大的计算机资源尚未在工程数值计算中获得应 用。因此我们使用分离式解法，v ，p 得代数方程的分离式解法的关键是 如何求解压力场，或者在假定了一个压力场后如何改进它。目前使用广泛 的s i m p l e 算法就是这么一种算法，s i m p l e ( s e m i i m p l i c i tm e t h o df o r p r e s s u r el i n k e de q u a t i o n s ) ，意即求解压力耦台方程的半隐方法，这一求解 不可压缩流场的方法是p a t a n k e r 与s p a l d i n g 在1 9 7 2 提出的。 一、s i m p l e 算法的基本思想： 在对n a v i e r - s t o k e s 方程的离散形式进行迭代求解的任一层次上，可以 给定一个压力场，它可以是假定的或是上一层次计算所得出的。一个正确 的压力场计算应该使计算得到的速度场满足连续性方程。但根据这样给定 的压力场计算而得的速度场，未必能满足连续性方程，因此要对给定的压 力场作改进，即进行修正，原则是：与改进后的压力场相对应的速度场能 满足这一迭代层次上的连续性方程。据此来导出压力的修正值与速度的修 正值，并以修正后的压力与速度开始下一层次的迭代计算。 s i m p l e 算法的计算步骤如下： 1 、假定个速度分布，记为u 0 , v o 以此计算动量离散方程中的系数及 常数项： 2 、假设一个压力场p ； 3 、依次求解两个动量方程得7 , ，v ； 4 、求解压力修砸值方程，得p 。； 5 、根据口改进速度值； 6 、利用改进后的速度场求解那些通过源项物性等与速度场耦合的 变量，如果西并不影响流场，则应在速度场收敛后再求解； 7 、利用改进后的速度场重新计算动量离散方程的系数，并利用改进 后的压力场作为下一层迭代计算的初值。重复上述步骤，直到获得收敛的 哈尔滨工程人学硕士学位论文 解。 二、s i m p l e 算法中引入了三个方面的假设： ( 1 ) 、速度场“o ，v o ，w o 的假设与p 的假设是各自独立进行的，两者间 无任何关系： ( 2 ) 、在导出速度修正值计算式时没有计及邻点速度修正值的影响； ( 3 ) 、采用线性化了的动量离散方程，即在一个层次的计算中，动量 离散方程中的各个系数及源项b 均假定为正值； 以上三点假设影响到速度场和压力场的协调和同步发展，因而会影响 收敛的速度；但是，由于可能存在的不协调性会随迭代的进行而逐渐减轻 以至消失，因而对最终结果并无影响。下面推导出压力和速度的校正方程。 1 、速度修正值的计算公式 为讨论方便，首先研究压力场改进后如何改进速度。设压力修正量， 速度修正量为p 。，“，v ，w 。，则改进后的压力与速度为： + p=p+p ： u = 1 , + “； v=v+v； 把p ，代入离散方程式得 a e ( ”。+ “，1 ) = 口柚( “柚+ + “n b 。) + 6 + 【( p p + p p 。) 一( p e + 圪) 】爿， “+ 是根据p + 值从上面方程解出的，满足 n 。“。= 口柚“。6 + 6 + ( p p 一p e * ) 4 。 这里我们假定由源项构成的b 的值不变，于是由上面两式相减得： l i e , 。= a n “。6 + ( p p 。一础4 假设在上式右边的两项中，压力修正值的影响是主要的，四周邻点的 速度修正值的影响可近似为零，于是得速度修正方程： l i e u 。= ( p p p f ) a 。 也可以写作为： “。= ( 争) ( p ；一p ：) = d 。( p ，一p 。) “p 于是改进后的速度为： “。= “。+ d e ( p p p e ) 哈尔滨工程大学硕士学位论文 类似的可得： v 。= v 。+ d 。( p p p ) 2 、压力修正值的代数方程 根据p 改进的速度场能满足连续方程，对于连续方程： 粤+ 掣+ 掣：0 ( 3 - 2 1 a ) a t敏a v 。 在控制容积上积分，有： 垃等业+ ( 川。一( 碱 6 y + ( 川nm ( 眺l a x = o ( 3 - 2 l b ) 整理后可的压力修正方程， n p p p2 d e p e + a w p - i - a s p s + n v p n + b t 3 - 2 1 c ) 其修正系数 a e = p 囊。6 y a w = p w d 。y a n = p n d n x d s = p s ds a x n p2 n # + n w + a n + d s 6 ：垡i 丝垒坐+ 【( 刖) 。一) 。+ 【( ) 。一) 。 (3a。i ( p ul ay(pvl a x 2 1 d ) 式( 3 - 2 1 e ) 就是确定压力修正值的代数方程。关于这一方程及其求解 要作以下几点说明： ( a ) 、如果速度场的当前值“，v 能使( 3 - 2 1 d ) 的右端等于零，则说明 该速度场已满足连续性条件，迭代业已收敛。因而6 的数值代表了一个控 制容积不满足连续性的剩余质量的大小。可以用各控制容积的剩余质量的 绝对值最大值，作为速度场迭代是否收敛的一个判据或指标。一种常用的 方法是以各控制容积6 的绝对值最大值及各控制容积的6 的代数和作为判 据，当速度场迭代收敛时，这两个数值都应为小量。 ( b ) 、根据p 计算而得的“，v 能使“= “+ “，v = v + v 1 满足连续性方程， 于是这样的“，v 就作为本层次上速度场的解，并且用它去改进离散方程系 哈尔滨上程大学硕士学位论文 数，从而开始下一层次的迭代计算。 3 、压力修正值方程的边界条件 为了求解式( 3 - 2 1 c ) ，还必须对压力修正值的边界条件作出说明。在 一般工程流场计算中，所见的边界条件是边界上 的法向速度已知或边界上压力分布已知。当边界 压力已知时，显然边界上的p 。t 0 ；我们这里详 细讨论一下边界上法向速度为已知的情形。在图 3 7 中，控制容积p 的u e 为已知，则对此控制容 积列出连续性方程时可以直接把已知的u 代入，e 不必再用( + 以) ，或者相当于以已知的代替 n 1 r 熨愚 一 帚心 爷 e 0 j s “：而令“：；0 。我们知道正是由于未知的：才需 图3 7 边界主控制容积 要引入p j 及p ：，既然u 。( 或者说u ：) 为已知，就不必再引入p ：。这相当 于在所形成的p j 的代数方程中口。p ；0 ，也就是* 0 。由此可见，无论 是边界压力为已知还是法向速度为己知，都没有必要引入关于边界上压力 修正值的信息。在计算中，可令与边界相邻的主控制容积的p 方程相应的 影响系数为零。 3 4 6s i m p l e c 算法( 1 9 8 4 ) 在s i m p l e 算法中，为了求解的方便，略去了速度修正值方程的 a n 。“二及a n h v n 。项，从而把速度的修正完全归结为由压差项的直接作用 所致。当我们在略去f d 。“：。等项不计时，实际上犯了一个“不协调一致” 的错误。为了能略去y 吼 “：。等项而同时又能使方程基本协调，对s i m p l e 算法进行修正得出了s i m p l e c ( s e m i i m p l i c i tm e t h o df o rp r e s s u r el i n k e d e q u a t i o n sc o n s i s t e n t ) 算法，在文献【2 0 j 中对两种二维层流流动的计算结果表 明，s i m p l e c 算法的特性远优于s i m p l e 。s 1 m p l e c 与s i m p l e 的计算 步骤基本相同，只有两点微小的区别： ( a ) 、以4 ( 吒一口) 代替s i m p l e 中的4 哎 哈尔滨工程人学硕士学位论文 ( b ) 、在s i m p l e c 算法中，p 不再亚松驰，即取口。= l