《气固稀相流中颗粒沉积和团聚的分子动力学模拟》由会员分享，可在线阅读，更多相关《气固稀相流中颗粒沉积和团聚的分子动力学模拟（8页珍藏版）》请在金锄头文库上搜索。
1、中国工程热物理学会 第t - 2 届学术会议 燃烧学 编号：0 6 4 1 7 5 气固稀相流中颗粒沉积和团聚的分子动力学模拟 李水清1J S M a r s h a U 2A R a t n e r 2 姚强1 1 清华大学热能工程系热能动力工程和热科学重点实验室北京1 0 0 0 8 4 中国 。2 爱荷华大学机械与工业工程系水利科学与工程所爱荷华城I A 5 2 2 4 6 美国 ( 联系电话：0 1 0 - 6 2 7 8 2 1 0 8 ；E - m a i l ：l i s h u i q i n g t s i n g h u a e d u c a ) 摘要本文在分子动力学模拟
2、基础上采用了考虑黏附作用的J K R 理论取代传统H e r t 作为颗粒 间作用力，并吸取了D E M 中滑动摩擦、滚动摩擦和扭动摩擦的思想建立一种细颗粒碰撞和团聚过程的 多时间尺度粒子动力学模拟方法。成功地实现了颗粒在纤维上沉积过程的模拟，结果表明微米颗粒的 V a nd e rW a a l s 力量纲约是曳力2 5 倍，只有中心流线附近颗粒更易于在纤维上沉积，进而解释了纤维 对微米颗粒捕捉效率较低的原因 关键字分子动力学模拟离散元模拟颗粒沉积纤维过滤低雷诺数流 分类号：T K l6 0 3 5 7文章标识码： A 1 前言 一 煤粉悬浮燃烧、水煤浆输送、颗粒过滤和受热面颗粒沉积等稀相流
3、动系统中，颗粒 问碰撞机理或近距离颗粒间微观作用( 吸引或排斥) 占有主导性作用。同时在燃烧碳烟形 成过程中颗粒碰撞和团聚也主导了其成长过程，进而影响到辐射传热、火焰温度和污染 排放等1 2 I 。虽然已有研究中国内外学者也采用了L a g r a n g i a n 系模拟并引入了颗粒碰撞率 概念，但并没有真正地将颗粒间微观作用力纳入到牛顿第二定律方程中渖射，因此对一些 特殊场合如颗粒纤维过滤、受热面颗粒沉积等过程无法实现微观尺度动力学模拟0 过去1 5 年在气固稠相流动领域中随着T s u j i ( 1 9 9 2 ) 、Y u 和X u ( 1 9 9 7 ) 、R h o d e s
4、 ( 2 0 0 0 和H r e n y ae ta l ( 2 0 0 5 ) 等学者不断发展，对于毫米尺度的大颗粒体系，离散元模型( D E M ) 因成功地将H o o k 定律或H e r t z 接触理论引入到了颗粒运动的牛顿第二定律方程中，而逐 步超过双流体模型成为了国际研究热点I “ I 。同时在分子尺度领域，将V a nd e rW a a l s 分 子间力作用( 可简化为L c n n a r d J o n e s 势) 考虑到牛顿第二定律方程的分子动力学模拟( M D ) 也得到了广泛应用。上述两种方法尚未成功地拓展到细颗粒( 微米或纳米尺度) 模拟的因 素大体有二：
5、( 1 ) 对于细颗粒而言黏附力作用增强使得H e r t z 接触理论不再适用，但黏附 力又无法简单地考虑为L e n n a r d - J o n e s 分子作用势；( 2 ) 气固稀相流是一个涵盖流体时间、 颗粒空气动力学响应时间和颗粒碰撞时间的多尺度计算，而仅基于颗粒碰撞时间尺度微 观模拟的计算量和耗时量都过于庞大。 本文在基于考虑黏附作用的J K R 接触理论的M D 基础上，同时引入了D E M 中对滑 动摩擦、滚动摩擦和扭动摩擦的概念，以锅炉尾部烟气中颗粒在单纤维上沉积过程为例， 对细颗粒沉积和团聚的特性开展了可自动优化计算时间步长的多尺度模拟。 收稿日期： 基金项目：国家自
6、然科学基金项目( N 0 5 0 3 0 6 0 1 2 ) ；国家重点基础研究发展规划项目( 2 0 0 2 C B 2 1 1 6 0 0 ) 5 9 5 2 考虑颗粒团聚效应的分子动力学模拟方法 “ 1 庙tcq)l劬 图1 颗粒接触空间示意图图2 微米尺度颗粒之间的作用力模式 传统稠相流动的D E M 模拟方法是基于H e r t 作为颗粒间作用力的假设，并在分 子动力学( M D ) 上发展起来的，亦称之为G D ( g r a n u l a rd y n a m i c s ) ，但其未考虑颗粒间黏附 和团聚效应。本文采用G D 模拟则是采用J K R 理论( J o h n s
7、 o n ，K e n d a l l 和R o b e r t s ，1 9 7 1 ) 作 为颗粒间作用力并在M D 基础上发展起来，特别适用于细颗粒的模拟1 9 1 。 ， 如图1 ，任一颗粒的平动和转动方程可依据牛顿第二定律分别表示为： 肌皇= B + B + L ，i d 。2 = M ，+ M ( 1 ) 口l四 式中，、Q 、m 和J 分别是颗粒的平动速度、转动角速度、质量和转动惯量：疋、 最和E 分别是流体对颗粒作用力、布朗扩散力和V a nd e rW a a l s 黏附力；M F 是流体作 用力矩，而M 。则代表了V a nd e rW a a l s 黏附力矩和颗粒碰撞
8、力矩之和。 2 1 流体作用力及力矩的计算 。 流体对颗粒的作用力包括曳力、S a f f m a n 升力、M a g n u s 力、附加质量力、压力梯度 力、重力和布朗力等。无量纲分析表明对密度比较大的细颗粒( z = p ，p , 一如就有作用。上述三参数表示为： 足= 赫磊= 燕，嘞= ( 竽厂 式中，为V a nd e rW a a l s 表面势能，其确定后文将详细讨论。 法向耗散力可表示为： = 一V 足n ( 1 0 ) 其中颗粒相对速度= v j 一巧一G Q + - q ) n ，法向耗散系数= a ( m k u ) “2 。据M a r s h a l l ( 2 0
9、 0 6 ) 的研究，法向弹性系数为k = ( 4 3 ) Ea C t ) ：口表示为碰撞恢复系数P 的函数。 2 2 2 滑动摩擦和切向力模型 颗粒切向接触力可采用C u n d a l l 和S t a c k ( 1 9 7 9 ) 提出 s p r i n g d a s h p o t s l i d e r 模型1 3 1 ， 瓦一耐= 一( j V s ( f ) 蟛) t s r r v 。t s ， ( 1 1 ) 右边第一项是弹性项，第二项是粘性摩擦；切向滑移速度v 。= v 。一( v 尼n ) n ，切向单位 矢量为t 。= v 。I v 。l ；切向弹性系数可表示为
10、砖= S G 口( f ) ，切向耗散系数珊= 佩。颗粒 ， 切向上还存在着滑动摩擦，粗颗粒流动中己共识符合库仑屈服准则，而对于细颗粒， T h o r t o ne ta 1 ( 1 9 9 1 ) i i ；l JV a nd e rW a a l s 力导致了其符合A m o n t o n 准则0 4 1 ，即 5 9 7 - - 纷1 r + 2 F c l ( 1 2 ) ，为滑动摩擦系数。根据滑动摩擦切向力的取值为C = m i n ( F ，, 一耐，) 。 2 2 3 滚动摩擦和扭动摩擦对转动方程修正 D o m i n i k 和T i e l e n s ( 1 9 9
11、5 ) 在太空尘埃聚合模拟中意识到存在着种滚动摩擦阻力使 的颗粒转动角速度减缓，并完善了相应计算模型1 2 l ；M a r s h a l l ( 2 0 0 6 ) 认为当颗粒转动角速 度在法向方向n 上不相同时还会产生扭动摩擦，并完善了其阻力矩的计算方程1 1 2 l 。具体 模型计算公式和参数列于表l 。 表1 滚动摩擦和转动摩擦阻力矩的计算公式 3 结果和讨论 y 图3 纤维捕捉细颗粒实验的示意图图4 不同时刻单纤维捕捉细颗粒的实验 3 1 “模型纤维”小R e 粘性流场的数值计算 笔者曾详细报道了采用显微成像技术所开展的单纤维捕捉细颗粒的实验l l 副，图3 是 捕捉实验的示意图，
12、图4 是不同时刻纤维上颗粒沉积图像，可观察到在纤维迎风面可形 成长度1 0 d 。的颗粒絮，说明了微米颗粒的团聚效应。 纤维过滤是一个经典的小雷诺数S t o k e s 流问题，为了能方便地求解流场，同时也为 可以在气溶胶过滤方面形成系统的研究体系，国外学者提出“模型纤维”的概念| 1 们，如 图5 所示，模型纤维就是指的气流垂直通过一排纤维的情况，= Z ，d ，为间距直径比， 在Y = + W 2 处气流满足周期性边界条件，因此，在周期性边界条件的作用下，一组纤维 情况就可以简化成了一个纤维来计算或模拟，当 1 时则可以视为单纤维。对于符合 周期性边界条件的模型纤维，比较经典的分析屑方法
13、是M i y a g i ( 1 9 5 8 ) 分析解和W a n g ( 2 0 0 4 ) 分析解1 1 7 1 ，此外对于特殊的单纤维还存在F r e e c e l l 模型，根据边界条件处理不同，大致 5 9 8 分为H a p p e l ( 1 9 6 3 ) 分析解和K u w a b a r a ( 1 9 5 9 ) 分析解，但事实上F r e e c e l l 模型是针对稠相 过滤体系的模型，这种近似对于单纤维和模型纤维都存在一定的误差捕l 。因此本文采用 C F D 数值求解的方法，具体使用F l u e n t 平台求解流场。无量纲化后流体控制方程为： 审d =
14、o ；p d 寺a = 彩+ p v 专2 a 流场无量纲参数为 ，a 2 u 。；文。x 勺= x L ；多= p 叼；p = v l r U o = R e l ，2 ( 1 3 ) ( 1 4 ) 由于Z 方向上纤维流场没有变化，为节省时间，流场采用二维计算，而颗粒场采用 三维计算。流场计算采用非结构化三角网格，为边界条件设置更合理并得到更精确流场 解，本文流场计算区间要大于颗粒场的计算区间，设= 1 ；则流体计算区间 - 2 夕2 ( W = 4 ) ；- 1 0 曼1 0 ( L = 2 0 ) 。 一$ 一一 图5 模型纤维示意图图6 颗粒计算区间Y 剖面内流线( 正交网格) 由于
15、三角网格计算结果是基于三角元顶点的，为得到颗粒所在位置处的精确流场， 需要采用P a r t i c l e i n c e l l 的搜寻和定位算法，判断颗粒具体在哪个三角单元里面。此外， 采用P o i n t ( p a r t i c l e ) - i n c e l l 思想将三角网格转化为正交网格是一种更有效率的方法I l 引，限 于篇幅本文不作深入探讨。图6 给出的是转化为正交网格后流体通过纤维的流线分布， 对应R e m = O 1 4 5 、，= 2 0 I m 和G o = 0 2 6 m s ，和K u w a b a r a ( 1 9 5 9 ) 分析解相比，数值解
16、结 果显示在夕= 2 处流线更平缓，而分析解求得流线更趋向中心对称线。 3 2 三维计算区间和无量纲颗粒模型参数确定 如图7 所示，颗粒计算区间为- - 2 三2 、一2 争2 和- 4 主4 ，在X Y 平面内小 于流体区间，与流体采用三角网格不同，颗粒计算采用正交网格且旃= 够= 出= o 1 ，目 前可实现5 万个颗粒计算。当想得到可视化较好的结果时，可以根据需要减少z 方向宽 度，并沿z 方向只布置一层颗粒。颗粒初始化可采用跟踪流线测试的方法，也采用随机 方法，Y 、Z 方向颗粒则采用周期性边界条件1 1 2 I ，显然气流x 方向采用颗粒流入一流出边 界条件更符合实际情况，本文旨在介绍颗粒与纤维表面的微观作用机理和动力学机制， 简化起见采用了周期性边界条件，而流入流出条件将另文讨论。 由前文已知流场参
