（安全技术及工程专业论文）考虑热对流的颗粒沉降的直接数值模拟.pdf

中北大学学位论文 d i r e c tn u m e r i c a ls i m u l a t i o no ft h ep a r t i c l e ss e d i m e n t a t i o nw i t h t h e r m a lc o n v e c t i o n a nk a n g s u p e r v i s o r ：p r o f c h a n gj i a n z h o n g a b s t r a c t o w i n gt ot h ew i d ee x i s t e n c eo ft h em u l t i p h a s ef l o w si nn a t u r a la n di n d u s t r i a lp r o c e s s e s ， t h er e s e a r c ho nt h ep h e n o m e n o n ，m e c h a n i s ma n dp r o c e s so ft h em u l t i p h a s ef l o w si so f i m p o r t a n ts i g n i f i c a n c e t h ea r b i t r a r yl a g r a n g i a n - e u l e r i a n ( a l e ) t e c h n i q u ew a su s e di nt h e d i r e c tn u m e r i c a ls i m u l a t i o no ft h ep a r t i c l e ss e d i m e n t a t i o nw i t ht h e r m a lc o n v e c t i o ni nt h e p a r a l l e lw a l l s t h ef l u i dm o t i o ni sc o m p u t e df r o mt h en a v i e r - s t o k e su s i n gt h ef i n i t e e l e m e n t m e t h o d t h ep a r t i c l ew a st r a c k e da c c o r d i n gt ot h ee q u a t i o n so fm o t i o no far i g i db o d yu n d e r t h ea c t i o no fg r a v i t ya n dh y d r o d y n a m i cf o r c e sa r i s i n gf r o mt h em o t i o no ft h ef l u i d ，a v o i d i n g l o t so fe x p e r i e n c ea n ds u p p o s i t i o no ff o r c ee n f o r c e do np a r t i c l e s t h r e es i t u a t i o n st h a ta r e p a r t i c l es e d i m e n t a t i o nw i t h o u tt h e r m a lc o n v e c t i o n ，s e d i m e n t a t i o nw i t ht h e r m a lc o n v e c t i o ni n h o to rc o l df l u i da r es t u d i e du s i n gt h ed i r e c tn u m e r i c a ls i m u l a t i o nm e t h o d t h en u m e r i c a l r e s u l t ss h o wt h a tt h et h e r m a lc o n v e c t i o na n dw a l le f f e c tw h i c hc o n t r o lt h ep a r t i c l e sm o t i o n a n di n t e r a c t i o n as i n g l ei s o t h e r m a lp a r t i c l ee x p e r i e n c e sd i f f e r e n tr e g i m e so fm o t i o nw i t ht h e s e d i m e n t a t i o nr e y n o l d sn u m b e ri n c r e a s i n g ：s t e a d ym o t i o nw i t ha n dw i t h o u to v e r s h o o t ；w e a k ， s t r o n ga n di r r e g u l a ro s c i l l a t i o n s c o m p a r e dw i t hi s o t h e r m a lp a r t i c l e ，t h et h e r m a lc o n v e c t i o n v a r i e st h ep a r t i c l es e d i m e n t a t i o nv e l o c i t ya n di t so s c i l l a t i o na m p l i t u d e ，a n dt h ep a r t i c l e m i g r a t e so f f t h ec e n t e r l i n ea tl o wr e d u r i n gt h es e d i m e n t a t i o no fc o l dp a r t i c l ei nt h eh o tf l u i d ， t h e r ea r ev o r t e xs h e d d i n ga r i s i n gi n d u c e db yt h et h e r m a lc o n v e c t i o n ，t h es e d i m e n t a t i o n v e l o c i t ya n dt h eo s c i l l a t i o na m p l i t u d eo fp a r t i c l ea r ei n c r e a s e d ，t h ep e r i o d i cd r a f t i n g ，k i s s i n g a n dt u m b l i n g ( d k t ) s c e n a r i ow a sf o u n di nt h e s et w oc o l dp a r t i c l e ss e d i m e n t a t i o n ，a n dt h e c o l dp a r t i c l e st e n dt od i s p e r s e d u r i n gt h es e d i m e n t a t i o no fh o tp a r t i c l ei nt h ec o l df l u i d ，t h e w a r mw a k ef o r m sas t r o n gu p w a r dt h e r m a lp l u m e ，t h es e d i m e n t a t i o nv e l o c i t ya n dt h e o s c i l l a t i o na m p l i t u d eo fp a r t i c l ea r ed e c r e a s e d ，t h ed r a f t i n g ，k i s s i n ga n dt u m b l i n gs c e n a r i o w a sn o tf o u n di nt h es e d i m e n t a t i o no ft w oh o tp a r t i c l e s ，a n dh o tp a r t i c l e st e n dt oa g g r e g a t e u n d e rt h es a m eh e a tt r a n s f e rc o n d i t i o n t h ei n f l u e n c eo ft h e r m a lc o n v e c t i o no np a r t i c l e 中北大学学位论文 m o t i o nw i l lb ed e c r e a s e dw i t ht h ep a r t i c l e f l u i dd e n s i t yr a t i oi n c r e a s i n g t h es u p p r e s s i o no f t h es e d i m e n t a t i o nv e l o c i t ya n dt h eo s c i l l a t i o na m p l i t u d eo fp a r t i c l eb yw a l le f f e c tw e r e i n c r e a s e dw i t ht h ec h a n n e lw i d t hd e c r e a s i n g k e yw o r d s ：a r b i t r a r yl a g r a n g i a n e u l e r i a n ，t w o p h a s ef l o w s ，d i r e c tn u m e r i c a ls i m u l a t i o n ， t h e r m a lc o n v e c t i o n ，s e d i m e n t a t i o n 原创性声明 本人郑重声明：所呈交的学位论文，是本人在指导教师的指导下，独 立进行研究所取得的成果。除文中已经注明引用的内容外，本论文不包含 其他个人或集体已经发表或撰写过的科研成果。对本文的研究作出重要贡 献的个人和集体，均已在文中以明确方式标明。本声明的法律责任由本人 承担。 论文作者签名： 安庭 日期： 2 归了奠z ，厂 关于学位论文使用权的说明 本人完全了解中北大学有关保管、使用学位论文的规定，其中包括： 学校有权保管、并向有关部门送交学位论文的原件与复印件；学校可 以采用影印、缩印或其它复制手段复制并保存学位论文；学校可允许学 位论文被查阅或借阅；学校可以学术交流为目的，复制赠送和交换学位 论文；学校可以公布学位论文的全部或部分内容( 保密学位论文在解密 后遵守此规定) 。 签 名：日期：力叨猡。f 夏广 导师签名：了。氓日期：沙固六修 中北大学学位论文 1 1 研究背景及意义 1 绪论 颗粒两相流动大量存在于自然界和工业设备中，如风沙泥沙输运、污水处理、防磨 润滑和能源动力、石油化工等工程领域中。对其流动现象、机理和过程的研究无论是在 学术理论还是在工程应用上都有重要意义，在其研究方法上即使一个微小的改进也可能 会在工程实践中产生巨大的经济效益。传统上，研究多相流体动力学的方法主要有两类： 一是分析方法，通过建立数学模型，经过大量的简化求得问题的分析解；二是实验方法， 通过实验，得到大量的实验数据，总结出半经验的公式来解决问题。在这些方法中，对 一些流动现象，如对颗粒的受力和颗粒的热量质量传递的计算，大多是在单个颗粒或固 定颗粒的基础上加以各种修正而得到的，许多实验也都是摸索和修正这些经验关系式， 不能刻画流场内部流动的细节。从分子运动论出发的微观分析法，由于物理上和数学上 的许多困难，目前还不能使用分子运动论来处理两相流动的实际问题。过去的工程设计 也主要依靠于经验公式，如颗粒的阻力系数( c d ) 、传热的努谢尔数( n u s s e l tn u m b e r ) 、 传质的沙伍德数( s h a w o o dn u m b e r ) 等的经验修正。由于这些经过修正的式子不能完全反 映出多种变化条件下的颗粒两相流动的本质特征，自然就会在工程设计中带来较大的误 差和浪费。随着计算机科学的迅猛发展，数值模拟越来越成为多相流研究的重要工具。 截止目前，虽然用于多相流的数值模拟方法有了相当大的进展，但在许多方面工作 仍然是初步的，例如，在数学模型中怎样考虑颗粒相和瞬态流动结构之间复杂的相互作 用、相变、化学反应等。现有的许多数值模拟方法，如均相模型、漂移模型、欧拉两流 体模型等由于采用了一些简化和假设，也很难捕获多相流场内部结构变化的细节和相问 相互作用，也无法深刻研究多相流动中存在的的各种相互影响和相互作用机制。随着人 们对多相流动现象认识的不断加深以及计算流体力学的飞速发展，使我们有可能对复杂 的多相流过程进行比较精确的数值模拟。这方面的工作既有重要的工业应用价值，又有 很高的学术意义，因而越来越受到人们的关注。 颗粒两相流动的微观流动现象和规律对于准确建立其宏观流动模型有着重要的意 1 中北大学学位论文 义。虽然直接模拟的计算量相当大，但它是获得多柑流动的微观流动现象与规律的有效 方法。以颗粒在流场中的受力为例，当采用模型化方法计算时，如欧拉一欧拉模型、欧 拉一拉格朗日模型等，在一些以特定的实验为基础或在一些假设后，可对某些问题给出 较满意的结果，但对大部分问题的计算这些模型都受到普遍性的限制，并不总是能得到 满意的结果，也很难捕获和深刻研究多相流场内部结构变化的细节和相间相互作用的机 制。产生该问题的主要原因之一是在多相流动问题中颗粒相的受力关系和其应力封闭还 没有一个完善的理论。以往在数值计算中颗粒所受到的气动力主要采用单颗粒的理论分 析或一些实验测量所得的结果，或在此基数上进行一些修正，如源于经典流体理论与实 验的单颗粒曳力关系在许多计算中被广泛采用。但是对于颗粒多相流问题，颗粒与流体 的耦合作用，颗粒间的相互影响等问题，必然使颗粒的受力变得复杂化，源于理想条件 下所获得的受力关系必然与所研究的问题存在一定的出入。另一方面，采用实验测量进 行有相对运动的颗粒间的受力研究是很难进行的，以往的颗粒受力实验都集中于固定颗 粒情况。而直接数值模拟却提供了一个有效的方法，它将计算网格缩4 , n 颗粒尺寸以下 进行流动的计算，颗粒占据有限体积，颗粒的受力不通过模型计算，而是通过积分表面 的粘性力和压力获得，所以称其为真正直接数值模拟( t d n s ) 。它可以给出实际流动中 颗粒的受力规律和颗粒的形状变化，也可以给出颗粒的脉动及其尾涡对流体的作用和颗 粒和流体之间的热量和质量的传递。因此，对多相流动进行直接数值模拟，所得到的研 究结果对于建立与实际情况相符的颗粒受力关系和揭示两相流动的相间相互作用的流 场特征和细节有非常重要的理论意义。它与以往的方法相比，为我们更准切地研究多相 流复杂流动现象开辟了新的道路。 两相流动力学在2 0 世纪6 0 年代开始迅速发展，并在8 0 一9 0 年代逐渐形成一个新 的学科，无论是研究方法和研究领域都有较大的发展。但如上所述，由于研究方法的限 制，人们对有传热传质化学反应和相变的多相流动的机理还未有很好的掌握和理解，其 困难在于基于颗粒相尺寸长度上的一些流场的介观特征，如在相间界面上的相互作用而 引起的流体流动和传热传质的耦合、相间界面的移动和变形以及颗粒相周围发生的反应 等，在很大程度上还未被充分的研究。这些问题和困难是多相流研究领域近期面临的挑 战。 为了更深入地研究颗粒两相流动的本质，从更小的介观尺度上( 如分散相单元的颗 2 中北大学学位论文 粒和液滴的大小) 阐明颗粒两相流动系统的本质特征和多相流场结构特征变化的机制， 本文将在a l e ( a r b i t r a r yl a g r a n g i a n e u l e r i a n ) 算法模拟等温颗粒两相流的基础上， 发展在热对流情况下对颗粒两相流进行真正直接数值模拟( t d n s ) 的方法，研究颗粒两 相流的运动规律和热对流对颗粒两相流运动的影响作用。 1 2 研究现状 对单个颗粒传热传质的早期研究，建立在外部流场己知且不受相间传送过程影响的 假定之上。a c r i v o s 和t a y l o r n l 在低p e c l e t 数下推导出绕球s t o k e s 流中的传热传质的 对称解，大p e c l e t 数下的解由l e v i c h 用边界层近似而获得瞳1 ，研究者通过数值求解的 方法也得到了任意p e e l e r 数下单个颗粒传热传质的情况口1 。 在单个颗粒的流动和传热传质的耦合方面，人们已经在液滴燃烧上做了一些工作 h 6 1 。这些研究中大都选用了薄火焰模型假定，燃烧反应被当作快速均一的化学反应。 大多数这样的研究中忽略了相间界面形状的变化，液滴被认为一直保持球状。 对复杂对流和分界面形状改变的研究集中在几个实验之中，这些实验都是研究固体 颗粒在热流体中的溶化和分解行为盯，引，多数研究者采用一些与溶解层有关的实验参数来 预测溶解的速度。实际上在此过程中，已观察到固体颗粒形成非均匀的表面形状，并且 发现相间界面形状对颗粒的传热传质有较大的影响。但很少有人用数值方法来模拟这一 过程。目前已得到的成果基本集中在对单个颗粒行为的研究上，它们对多相流动的自然 本质的理解没有太多的贡献。到目前仅有为数不多的几篇文章对颗粒移动和分界面变形 方面作了理论研究田j m l i , 1 2 。在这些研究中有的忽略了传热传质对流，有的只关心纯粹的 物质相变，没有考虑质量传递。另一方面，采用实验测量进行有相对运动的颗粒间的受 力和传热传质特性研究是很难进行的，已有的一些实验工作集中在固定颗粒的情况，而 数值模拟却提供了一个有效的方法。国内许多单位和研究人员，如清华大学、中科院、 西安交通大学、浙江大学等对颗粒两相流动的数值模拟作了许多卓有成效的工作，取得 了许多重要成果n 3 j4 1 5 t 1 6 1 7 t 1 8 旧2 0 2 “2 2 2 3 卫4 _ 5 矧。但在数值模拟方面，如何应用合适的方法在 数学模型中考虑颗粒相和瞬态流动结构之间复杂的相互作用以及传热传质相变、化学反 应等方面还需进一步探索和研究。 对两相和多相流动进行完整的数值模拟所需知识的深度和广度是惊人的，在不懈的 3 中北大学学位论文 探索过程中人们提出了多种数理模型，从最简单的均相流模型、分相流动模型、漂移模 型，一直到双流体模型。以k 模型中，均相模型完全没有考虑两相的差异，分相流动模 型和漂移模型在一定程度上引入了两相的相互作用，但仍过于简单而无法精确描述两相 的运动与空间分布。在目自订应用的双流体模型中，针对两相分别写出质量、动量和能量 守恒方程，在基于实验现象为基础的这类模型中流体和颗粒问的质量、动量和热量的传 递都须用经验公式输入，对特定的一些与实验相似的情况，如一些循环流化床中的多相 流动的颗粒运动、传热、传质进行了较好的数值模拟，但在模型的应用上受到了普遍性 的限制。在基于流体力学的这类模型中，由于颗粒的粘性应力封闭还没有一个完善的理 论，对相间起较大作用的曳力项，由于难以引入反映颗粒群结构的参数，其计算时的系 数都是在单颗粒曳力系数基础上修正得到的，而且做了许多理想化的假设，忽略了许多 应该考虑的因素心7 施巩弛3 “3 2 。引。在具体的数值模拟中大都采用了欧拉两流体方法，丢失 了许多真实的离散相流动图象，即使采用欧拉一拉格朗日方法，通常对颗粒的计算也采 用了一些经验系数，源于经典流体理论与实验的单颗粒曳力关系也在许多计算中广泛使 用，同样避免不了要丢失一些真实的离散相流动图象。但是对于颗粒两相流动问题，颗 粒和流体的紧密祸合作用、颗粒间相互影响等问题，必然使颗粒的受力变得复杂化，源 于理想条件下的颗粒受力关系必然与所研究的实际问题存在出入，影响数值模拟的效 果。因此这些方法大多不能精确描述离散相和流动结构之间复杂的相互影响机制，流动 的介观特征没有被真正考虑在内。流体拟颗粒模型和无网格粒子法( s p h ) 理论上可以 解决其它模型所无力解决的微观流场结构问题。然而，将宏观连续的介质离散化带来的 复杂性是巨大的，限于计算能力和对微观运动方式的理解，短期内在实际当中应用还有 较大困难啪伪3 0 3 1 3 2 ? 。截止目前，由于真正直接数值模拟热对流情况下颗粒两相流的复 杂性，国内外在这方面的研究，特别是国内用t d n s 的方法模拟相变热对流条件下的颗 粒两相流发表的文献还较少见到。 综上所述，对单个颗粒和固定颗粒的两相流动、传热传质已作了大量的研究，但对 热对流情况下颗粒两相流动的本质特征到目前为止很少被涉及。在目前许多的数值模拟 中，采用的许多模型和方法都做了大量各种简化，无法得到真实的离散相流动图象，忽 略了离散相的形状变化及其在相变热对流条件下对连续相的反作用，不能精确描述颗粒 和流动结构之间复杂的相互影响机制，流动的介观特征没有被真正考虑在内。本文将在 4 中北大学学位论文 对等温颗粒两相流进行直接数值模拟的基础上，发展在热对流情况下对颗粒两相流进行 直接数值模拟的方法，应用该方法，在牛顿流体中跟踪颗粒的运动和表面形状的改变， 强调分界面上相间相互作用和流动的介观特征，捕获因相间的相互作用引起的多相流动 的特征变化等流场内部结构的细节，获得并归纳总结其动力学和传热学特性，深入研究 有热对流条件下颗粒两相流的流动本质。 1 3 研究方法 本文应用任意拉格朗日一欧拉( a l e ) 方法对颗粒两相流进行真正直接数值模拟， 采用有限元方法数值求解流场的n - s 方程，并增加联立求解能量方程，应用牛顿定律跟 踪颗粒运动，并通过积分颗粒表面的粘性应力和压力获得颗粒的受力，避免了其它模型 中对颗粒受力采用的许多假设，从而实现了对颗粒两相流运动的真正直接数值模拟。数 值模拟的网格缩小到颗粒尺寸以下，颗粒占据有限体积，通过d e l a u n a y - - v o r o n o i 法生 成非结构化的三角形单元网格；颗粒移动时，通过求解l a p l a c e 方程得到网格移动速度， 当单元网格严重变形时网格将重新划分，以确保网格质量。颗粒受力方程、扭矩方程和 流体的动量方程一起进行有限元的g a l e r k i n 法推导，这样颗粒和流体间相互作用的力 和扭矩就不必专门加以计算。颗粒位置的更新将由其速度决定，时间步长由颗粒的速度 和加速度来自动调整，方程的非线性部分由牛顿迭代求解，线性部分由g m r e s 算法来求 解。 本文对所研究对象的数学描述上，对连续相来说其控制方程为质量、动量、能量守 恒方程，密度和粘性的变化考虑流场温度变化的影响，热对流通过自然对流、强制对流 和混合对流来考虑，颗粒相的位置和形状及界面温度为连续相的边界条件，连续相在分 界面的速度即颗粒相表面各点的速度，颗粒的运动由牛顿第二定律控制。相间耦合通过 颗粒和流体间相互的作用力和力矩及边界条件的施加实现。 1 4 本文的研究内容 本文在前人的基础上开发并完善了模拟程序，使用a l e 算法对考虑热对流的颗粒在 通道中沉降的过程进行真正直接数值模拟( t d n s ) 。模拟计算在牛顿流体中的颗粒两相 5 中北大学学位论文 流动，得到其介观尺寸上的流场结构，捕获因相问的相互作用引起的多相流动的特征变 化等流场内部结构的细节，获得并归纳总结其动力学和传热学特性，与相关实验和数值 模拟结果进行比较，并阐明其流场内部结构及相间的耦合与相互作用的机制。具体如下： ( 1 ) 分别模拟了在通道中绕圆柱的强制对流和同心圆环之间的自然对流，做为测 试算例。 ( 2 ) 模拟了考虑热对流的单颗粒在竖直通道中的沉降过程，分别模拟了颗粒在等 温流体、热流体和冷流体中的沉降，并进行分析。 ( 3 ) 模拟了不同颗粒和流体密度比情况下的三种沉降过程，进行对比，并分析。 ( 4 ) 模拟了不同通道壁宽度情况下单颗粒的沉降过程，并分析。 ( 5 ) 模拟了双颗粒分别在等温流体、热流体和冷流体中的沉降过程，并分析。 6 中北大学学位论文 2 直接数值模拟理论及控制方程 数值模拟方法在两相流研究中越来越具有重要的地位，人们提出了多种模型和格式 来数值模拟两相流流场，但对两相流流场的直接数值模拟只是最近时期随着计算机技术 的发展才逐渐开展。直接数值模拟( d i r e c tn u m e r i c a ls i m u l a t i o n ，d n s ) 是一种不引 入模型的数值计算过程1 。 2 1 直接数值模拟方法概述 所谓两相流的直接数值模拟就是把流场中的颗粒周围计算网格缩小到颗粒尺寸以 下进行流动的计算，颗粒占据有限体积，颗粒的受力不通过模型进行计算，而是通过积 分表面的黏性力与压力获得。所以此方法又称为真正直接数值模拟( t r u ed n s ，t d n s ) 。 如图2 1 所示，直接数值模拟中使用的非结构化网格，其中颗粒周围的网格尺寸要小于 颗粒尺寸。 图2 1 直接数值模拟中的非结构化网格 要进行两相流的d n s 计算，必须同时求解控制流体运动的n a v i e r s t o k e s ( n - s ) 方 程以及控制颗粒运动的刚体运动方程。这两类方程之间是通过颗粒边界的无滑移条件以 及出现在刚体运动方程中的流体作用力和力矩耦合起来的。固体颗粒的受力通过积分流 体作用在颗粒表面的黏性力和压力获得，而不是通过模型来计算。 在两相流中，牛顿流体的n s 控制方程如下： v ，一0 ( 2 1 ) 7 中北大学学位论文 p 。皆矿w ) 印川2 昭亿2 ， 假设颗粒为球形，颗粒质心的坐标为x ( f ) ，颗粒运动方程如下： ，l f 盟：g f + r d t ( 2 3 ) ，f 丝。乃 d t ( 2 4 ) 一d x 。m ( 2 5 ) d t 、。 塑；q ：( 2 6 ) i 表示第i 个颗粒，t 表示时间，m i 和 分别是颗粒的质量和转动惯量，m 和f 2 i 是颗 粒的线速度和角速度，g f 、r 和a 分别是体积力，流体对颗粒的作用力和力矩。其中互 是流场对颗粒的作用力，互是流场对颗粒的作用力矩。x 和p 分别表示颗粒位移和旋转 角度。在固体颗粒表面，流体的速度满足固壁无滑移条件： ，= 杉+ q fa ， ( 2 7 ) 作用在颗粒表面的力e 是颗粒表面所受拽力的总和： r = 盯万一一p 托+ 2 皿，l ( 2 8 ) 其中b ：要( v y + v v r ) 是流体变形率，玎是外法线单位向量。 2 2 任意拉格朗日一欧拉方法 非线性连续介质力学的有限元分析方法己取得了很多令人瞩目的进展，其中在大变 形问题中绝大多数研究工作都采用拉格朗日( 主要用在固体力学中) 或欧拉( 主要用在 流体力学中) 描述方法口引。 在拉格朗日描述( 即物质描述或l 描述) 中，计算网格固定在物体上随物体一起运 动，即网格点与物质点在物体的变形过程中始终保持重合，因此物质点与网格点之间不 存在相对运动( 即迁移运动，也称对流运动) 。这大大地简化了控制方程的求解过程， 8 中北大学学位论文 而且能准确描述物体的移动界面，并可跟踪质点的运动轨迹。但在涉及到特大变形的问 题中，物质的扭曲将导致计算网格的畸形而可能使得计算失败。 在欧拉描述( 即空间描述或e 描述) 中，网格固定在空间中，即计算网格在物体的 变形过程中保持不变，因此可很容易处理物质的扭曲。但对运动界面需要引入非常复杂 的数学映射，将可能导致较大的误差。另外当使用普通的伽辽金离散时，由于迁移项的 影响，有限元方程中的系数矩阵是非对称的，而且还可能得到振荡解。 纯拉格朗日和纯欧拉描述都存在严重的缺陷，但也具有各自的优势。如果能将二者 有机地结合，充分吸收各自的优势，克服各自的缺点，则可解决只用纯拉格朗日和纯欧 拉描述所解决不了的问题。任意拉格朗日一欧拉( a r b i t r a r yl a g r a n g i a n e u l e r i a n ，a l e ) 方法就是基于此目的最早由n o h ( 1 9 6 4 ) 以耦合欧拉一拉格朗日的术语提出的，并用有 限差分法求解带有移动边界的二维流体动力学问题。在n o h 的研究工作中，网格点可以 随物质点一起运动，但也可以在空间中固定不动，甚至网格点可以在一个方向上固定， 而在另一个方向上随物体一起运动，因此a l e 描述也被称为耦合欧拉一拉格朗日描述。 例如在液体表面波的传播问题中，网格点在垂向随物质点一起运动，而在水平方向上固 定不动。这样可很容易描述液体表面的运动，而且网格不会发生扭曲。 在a l e 描述中计算网格可以在空间中以任意的形式运动，即可以独立于物质坐标系 和空间坐标系运动。这样通过规定合适的网格运动形式可以准确地描述物体的移动界 面，并维持单元的合理形状。类似于e 描述，在a l e 描述下的控制方程中也将出现对流 项，因此也可能得到振荡解，需要进行相应的数值处理。纯拉格朗日和纯欧拉描述实际 上是a l e 描述的两个特例，即当网格的运动速度等于物体的运动速度时就退化为拉格朗 日描述，而当网格固定于空间不动时就退化为欧拉描述。 在流体一结构相互作用问题中，流体的运动幅度比较大，因此在固体域中一般使用 拉格朗日描述，而在流体域内使用欧拉描述。在流体一结构交界面上的每个点处需设置 两个节点，即流体节点和固体节点。整个变形过程中，在流体一结构交界面上流体质点 不能穿越网格，而且也不能穿越固体表面，同时应使固体网格点和流体网格点在界面上 始终保持重合。蝻t3 7 观3 。 9 中北大学学位论文 2 3 基于网格重构的d n s 方法 如图2 1 所示，颗粒在流体中的运动必然导致在每一个时间步内都需要进行网格的 重新划分，这种方法称为基于网格重构的d n s 方法b 3 1 。当然，在每两个时间步之间还需 要进行数据的投影以解决控制方程空间离散形式变化的问题。 为了解决重新划分网格带来的巨大工作量，在处理颗粒运动时需要使用任意拉格朗 日一欧拉( a l e ) 移动网格法h 仉4 “42 】，假设颗粒表面的网格结点跟随颗粒一起运动，而位 于流体中的网格结点位置则通过修正过的l a p l a c e 方程计算出来，以保证它们是光滑分 布的。在每一个时间步，网格结点均需要根据颗粒的运动来更新，一旦更新过的网格中 有的单元过于扭曲，那么就需要对整个流场中的网格重新划分。a l e 方法能够节省计算 时间，采用a l e 法需要将流场中的参数投影到新的网格系统上，颗粒位置、运动方向以 及网格点的位置都被显式地更新，而流体速度、颗粒的移动速度和角速度则被隐式地更 新。 a l e 方法在发展过程中又分为2 种类型，一种是整体方案( i n t e g r a t e dm e t h o d ) h 3 1 ， 一种是分裂方案( s p l i t t i n gm e t h o d ) h 4 1 。在整体a l e 方法种，速度场和压力场是同时 求解的，而在分裂a l e 方法种，速度场和压力场是分开按顺序求解的。与整体a l e 方法 相比，分裂a l e 法可以得到一个更小的方程组，但是，分裂a l e 法并不是在每一个时间 步长都将无散度条件强加入方程中，从而导致了速度场并不精确的满足动量方程，因此 必须使用比较小的时间步长来求解动量方程以减小误差。本文对颗粒两相流进行直接数 值模拟采用了整体a l e 方法。 2 4 颗粒沉降过程直接数值模拟的控制方程 在有牛顿流体的竖直通道中，释放一个初始速度为零的固体颗粒，由于颗粒密度略 大于流体密度，颗粒受到重力作用丌始沉降n 5 4 6 4 7 48 4 9 3 。假设颗粒是圆形的，流体初始温 度为丁o ，颗粒温度为l ，沉降过程中颗粒温度恒定不变( 这是由于我们研究的是热对 流对颗粒和流场的影响作用，如果颗粒温度不固定，很快就会与流体温度一样，热效应 的影响很快就会消失) 。如图2 2 所示，x 轴方向竖直向下，y 轴方向水平向右，颗粒的 初始位置为【u ，y o ，。 l o 中北大学学位论文 图2 2 计算区域不恧图 流体的连续性方程、运动方程、能量方程分别为： v v = 0 ( 2 9 ) p 。驴w ) 一唧川2 昭弦 叫。w p z亿 流体密度p 的变化与温度的关系为： p ；p o 1 一仃一t o ) 】 ( 2 1 2 ) 其中p o 是z o 温度时的流体密度，卢是流体的热膨胀系数。 颗粒运动方程如下： 肌f 坐。g i + r d t ( 2 1 3 ) f i d q _ _ _ 兰。死 d t ( 2 1 4 ) i 表示第i 个颗粒，m i 和 分别是颗粒的质量和转动惯量，1 , 5 和q i 是颗粒的线速度 和角速度，g i 、r 和n 分别是体积力，流体对颗粒的作用力和力矩。 中北大学学位论文 计算区域的四周的边界条件如下： b c ，c d ，d a 边j 二：v = 0 ，t = t o ： a b 边 州刎，芸- i - 。 f 表示切应力，在模拟过程中，计算区域将随着颗粒的运动而移动，并在每一时间 步内重新定义，但区域的大小和颗粒在区域中的相对位置不变，即颗粒在计算区域中与 顶端和底端的距离始终不变。在本文中所，和五分别为印z 2 4 和印矗4 3 2 。 采用如下无量纲参数： 尺e ：p o u d 雷诺数： a ( 2 1 5 ) g r 。p ；, p a ：r _ d 一 g 格拉晓夫数： ( 2 1 6 ) p r 。丝 普朗特数： k ( 2 1 7 ) 瑞利数：r a g 脚( 2 1 8 ) 在处理多颗粒之间、颗粒与通道壁之间的碰撞问题时，本文采用了如下方法：在颗 粒周围很小的范围内定义了一个安全区域，当颗粒间距离小于安全距离时，阻止颗粒接 近的排斥力被激活，排斥力将颗粒从安全区域推出，从而使颗粒与其他颗粒或与通道壁 分开。 假设& 和励为两圆形颗粒，半径分别为r i 、碍，质心分别为g i 、g ! ，则两颗粒b i 和马之间的排斥力气满足： ， f 略 置+ b + p f 2 1 三 ( g f g ，) 【( r + r j + p ) 2 - d ；】，i f d 。s r + 尺，+ p 2 1 9 其中：d q2 iq g jl ，p 是产生排斥力的范围，为“刚度”参数。 颗粒与通道壁之间的排斥力譬满足： 1 2 中北大学学位论文 f o ，i f d l 2 r + j d 髟2 土( g j 一咖2 r 训z 一珊m 啦+ p q 2 1 其中：d t ；lg f ql 是颗粒& 的中心与通道壁r ，另一侧假设颗粒的中心之间的距离 ( 如图2 3 ) ，为另外的刚度”参数。 f 。 i i d i r 、 ， l 、 ：7 图2 3 计算颗粒与通道壁之间排斥力时的假设颗粒 本章主要对直接数值模拟的理论进行了概述，介绍了真正直接数值模拟( t d n s ) 和 任意拉格朗日一欧拉( a l e ) 方法，并描述了本文编写程序所采用的两相流模型的n s 控制方程。 1 3 中北大学学位论文 3 测试算例 颗粒与流体有温差的两相流体中，两相流流场的运动要受到自然对流和强制对流的 影响，因此在数值模拟颗粒两相流之前，为了验证本程序对所模拟问题的准确性和有效 性，我们选择了两个测试算例，即在通道中绕圆柱的强制对流和同心圆环之间的自然对 流。两算例都已经被充分的研究，并有大量相关的参考文献，对本程序的模拟有重要的 参考意义。前者能够验证程序对强制对流模拟的准确性，后者则能验证程序处理自然对 流的能力。 3 1 通道中绕圆柱体的强制对流 c h a n g f i n l a y s o n ( 1 9 8 7 ) 和d e n n i s ( 1 9 8 6 ) 咖5 13 模拟了通道中绕圆柱体的强制 对流，本文使用相同的几何参数对该算例进行了模拟。 图3 1 、图3 2 分别为雷诺数2 0 时，圆柱体表面的努塞尔数( n u ) 分布和圆柱体表 面的压力系数分布，其中： 肌。丝。翌：里! 丝 ka t ； 。= 击 本文的方法计算得到的结果与c h a n g & f i n l a y s o n ( 1 9 8 7 ) 和d e n n i s ( 1 9 8 6 ) 计算 的结果相当的一致。 1 4 中北大学学位论文 图3 1 圆柱体表面的努塞尔数分布( r e = 2 0 ) 图3 2 圆柱体表面的压力系数分布 3 2 同心圆环之间的自然对流 为了验证程序处理自然对流的能力，本文采用了与k u e h n & o o l d s t e i n ( 1 9 7 6 ) 瞄刃对 此问题研究时相同的参数：设内圆环半径为忍，温度为乃，外圆环半径为忍，温度为l ， 15 中北大学学位论文 r o r i = 2 6 ，p r 一0 7 ，r a = 5x 1 0 4 ，其特征长度为m 一尼。 图3 3 、图3 4 、图3 5 分别为沿半径剖面上无量纲的温度分布、无量纲的速度分 布、同心圆环间的等温线和流线。其中： t ， l ，宰( r o r i ) y 一= 一 a a 口g 一 热扩散率：p c 可以看出本文的结果与k u e h n g o l d s t e i n ( 1 9 7 6 ) 的结果非常一致。 图3 3 沿半径无量纲温度分布( 0 = 0 、3 0 、9 0 、1 8 0 ) 1 6 中北大学学位论文 l 句 。 图3 4 沿半径无量纲速度分布( 口= 1 2 0 ) 10o o ， f ，“； 撼； 攥涝 德、 -暑碜一 l、n。，一 图3 5 同心圆环间的流线和等温线 17 中北大学学位论文 本章通过两个测试算例，模拟了有自然对流和强制对流的两相流流场，模拟结果与 相关参考文献相当的一致，证明了本文所采用的数值方法和所编写的程序可以证确模拟 考虑热对流的颗粒两相流的流场运动。 1 8 中北大学学位沦文 4 单颗粒在通道中的沉降 为了研究单颗粒在竖直通道中的沉降过程和热对流对单颗粒沉降的影响作用，本文 分别模拟了三种不同情况的颗粒沉降过程：颗粒和流体无温差的等温颗粒沉降、冷颗粒 在热流体中沉降和热流体在冷流体中沉降。 4 1 等温颗粒的沉降 在含有牛顿流体的的竖直通道中，将初始速度为0 、半径为d 的圆形颗粒由通道中 心线上释放，其中竖直通道宽l = 4 d ，颗粒中心与计算区域的顶端和底端的距离均为7 5 d 。 流体初始温度为r o ，颗粒温度为兀，设t o = 冗。颗粒的初始位置为( 0 ，y o ) ，颗粒密度p s 略大于流体密度p o ( 恒定不变) ，受到重力作用释放后开始沉降。通过改变颗粒的密度， 可控制颗粒x 方向达到的最终平衡速度，从而得到不同沉降情况的雷诺数间隔范围。 情况a ：( 0 1 r e 2 ) 图4 1 显示了0 1 r e 2 时，颗粒从不同水平位置释放后的运动轨迹。若颗粒在竖 直通道的中心被释放，颗粒将沿通道中心线稳定的沉降；当颗粒的初始位置不在通道的 中心时，颗粒由于通道壁的作用，释放后则向中心线方向运动，当颗粒到达中心线时， 便沿着中心线竖直稳定的沉降，并且雷诺数越大，颗粒达到中心线位置的时间越短。 1 9 中北大学学位论文 图4 1 情况a ：不同密度颗粒从不同水平位置释放后的运动轨迹。( a ) r e 2 ( b ) r e = 0 6 1 2 ( c ) r e = 1 2 2 情况b ：( 2 r e 2 1 ) 颗粒最终的平衡位置仍为竖直通道的中心线，但颗粒的运动轨迹不在是单调的曲 线，颗粒释放到达中心线后，将继续在水平方向上运动，然后形成了阻尼运动的轨迹， 最终稳定在中心线上( 如图4 2 ) 。 图4 2 情况b ：颗粒的运动轨迹。( a ) r e = 3 ，4 7 ( b ) r e = 1 5 8 2 0 中北大学学位论文 情况c ：( 2 1 r e 5 5 ) 当2 1 r e 5 5 时，颗粒若在通道中心位置释放，仍沿中心线保持稳定的沉降；当颗 粒不在中心线时，颗粒向中心线位置运动的同时，产生水平方向上的摆动，且摆动幅度 越来越小，最终沿着中心线稳定沉降( 如图4 3 ) 。 图4 3