（机械电子工程专业论文）基于格子boltzmann方法的海流能水轮机大雷诺数模拟研究.pdf

基于格子b o l t z m a n n 方法的海流能水轮机 大雷诺数模拟研究 捅晏 水轮机是海流能发电系统装置的重要组成部分，不同的水轮机的换能效果 也不相同，因此水轮机的设计优化对于提高海流能的利用效率具有至关重要的 作用。水轮机叶片数量的选择以及水轮机叶片翼型的水动力特性分析是水轮机 参数优化的基础工作。 格子b o l t z m a n n 方法是近几十年来发展起来的一种模拟复杂流体现象的介 观方法，近年来引起了国际上许多学者的关注，其理论不断成熟，应用领域不 断扩大，该方法相对于传统的数值模拟方法具有算法简单、精度高、可以直接 求解压力、能够模拟复杂边界条件、易于并行计算等优点，已逐渐发展成为仿 真模拟流体动力学问题的一类重要的数值方法。但是格子b o l t z m a n n 方法在模 拟宏观n s 方程时，为降低误差导致了该方法模拟大雷诺数流场的局限性，因 此大雷诺数下流场的模拟成为将该方法应用于实际工程领域的瓶颈。 本文将格子b o l t z m a n n 方法应用于海流能发电水轮机的水动力特性分析和 研究中，并对格子b o l t z m a n n 方法的基础性理论做进一步研究分析，通过基于 大涡模拟的格子b o l t z m a n n 方法( l b m l e s ) 解决了格子b o l t z m a n n 方法不能 模拟大雷诺数流场的问题，并通过数值模拟和对比分析，验证了基于大涡模拟 的格子b o l t z m a n n 方法模拟大雷诺数流场的准确性，从而为格子b o l t z m a n n 方 法应用于实际工程领域的海流能水轮机模拟奠定了基础；对于叶片翼型的研究 通过对不同攻角、不同厚度和弯度的翼型进行模拟对比分析，以及处理计算数 据总结了大雷诺数下不同翼型参数与翼型水动力特性的变化规律关系，并绘制 了相应流线图和水动力特性曲线等；另外还对整体的水轮机的绕流模拟进行了 模拟研究，模拟了不同叶片数量的垂直轴水轮机( 两叶片、三叶片、四叶片) ， 求解了不同叶片数量的水轮机的转矩系数与尖速比的关系，得出了不同叶片数 量的水轮机获能效率的一般规律，为水轮机的优化设计提供了理论依据。 关键词：格子b o l t z m a n n 方法；大涡模拟；叶片翼型；水轮机；水动力特性 r e s e a r c ho nl a t t i c eb o l t z m a n ns i m u l a t i o no nh i g hr e y n o l d s n u m b e rf l o wa r o u n dm a r i n ec u r r e n tt u r b i n e a b s t r a c t m a r i n ec u r r e n tt u r b i n eh a sa 啊t a li m p a c to nt i d a lc u r r e n te n e r g yc o n v e r s i o n d e v i c e t h ep o w e rc o e f f i c i e n to fm a r i n ec u r r e n tt u r b i n ev a r i e sw i mt h es t r u c t u r eo f t h et u r b i n e a st h eb a s i cr e s e a r c hw o r ko ft h eo p t i m i z a t i o no fm a r i n ec u r r e n t t u r b i n e s ，h y d r o d y n a m i c sa n a l y s i so fh y d r o f o i lb l a d e sa n dt h en u m b e ro fh y d r o f o i l b l a d e sa r ev e r yi m p o r t a n ti no r d e rt oi m p r o v et h ep o w e rc o e f f i c i e n to fm a r i n e c u r r e n tt u r b i n e s g r a d u a l l yd e v e l o p e di nr e c e n ty e a r s ，l a t t i c eb o l t z m a n mm e t h o dh a sd r a w n m u c ha t t e n t i o na n db e c o m ea l li m p o r t a n tm e s o s e a l es i m u l a t i o nm e t h o df o rc o m p l e x f l u i dm o d e l i n ga n ds i m u l a t i o n i th a s t h ea d v a n t a g e so fs i m p l ea l g o r i t h m ，h i g h a c c u r a c y , d i r e c t c a l c u l a t i o nf o rt h ep r e s s u r e , s u i t a b l ef o rc o m p l e xb o u n d a r y c o n d i t i o n sa n dp a r a l l e lc o m p u t a t i o n s b u ts t a n d a r dl a t t i c eb o l t z m a n nm e t h o di s u n a b l et os i m u l a t eo f h i g hr e y n o l d sn u m b e rd u et ot h ed e m a n do fe r r o r , w h i c hl e a d t ot h el i m i t a t i o n so f t h i sm e t h o do ne n g i n e e r i n ga p p l i c a t i o n i nt h i sp a p e r , l a t t i c eb o l t z r n a n nm e t h o di sa d o p t e dt ot h ea p p l i c a t i o nr e s e a r c h o nh y d r o d y n a m i c sa n a l y s i so fm a r i n ec u n e n tt u r b i n e t h r o u g hf u r t h e rr e s e a r c ho n b a s i ct h e o r yo fl a t t i c eb o l t z m a n nm e t h o d , l a r g ee d d ys i m u l a t i o nc o m b i n e dw i t h l a t t i c eb o l t z m a n nm e t h o di su s e dt os o l v et h ep r o b l e mo fs i m u l a t i o no fh i g h r e y n o l d sn u m b e r t h o u g hn u m e r i c a ls i m u l a t i o na n dc o m p a r a t i v ea n a l y s i s ，i ti s p r o v e dt h a tt h el b m l e si sa c c u r a c ya n de f f e c t i v eo ns i m u l a t i o no fh i g hr e y n o l d s n u m b e r b ys i m u l a t i n gf l o wa r o u n dh y d r o f o i lo nt u r b u l e n tu n d e rt h ec o n d i t i o n so f d i f f e r e n ta t t a c ka n g l e s ，d i t t e r e n th y d r o f o i lb l a d e sw i t hd i f f e r e n tc a m b e ra n d t h i c k n e s s , t h et r e n do f h y d r o d y n a m i c sp a r a m e t e f si nd i f f e r e n th y d r o f o i le o n d i t i o mi s s u m m a r i z e d a n dt h ec o r r e s p o n d i n gs t r e a m l i n e sa n dh y d r o d y n a m i cc u r v e sa g od r a w b yp r o c e s s i n ga n da n a l y z i n gt h ec a l c u l a t i o nd a t a b a s e do na b o v ew o r k , t h et o t a l t u r b i n ei sa l s os i m u l a t e d t h r o u g hs i m u l a t i o no fv e r t i c a lt u r b i n e s 丽t l ld i f f e r e n t b l a d e sn u m b e ra n dc a l c u l a t i o no fh y d r o d y n a m i c sp a r a m e t e r s ，t h er e l a t i o n sb e t w e e n n u m b e ro fb l a d e sa n dh y d r o d y n a m i c so ft u r b i n ea r ek n o w na n dr e g u l a rp a t t e r no f p o w e rc o e f f i c i e n tv a r i e sw i t hn u m b e ro fb l a d e si so b t a i n e d t h er e s e a r c hw o r ki n t h i sp a p e rp r o v i d e sr e f e r e n c ea n dg u i d a n c ef o rd e s i g n i n ga n do p t i m i z i n gm a r i n e k e y w o r d s ：l a t t i c eb o i t z m a n nm e t h o d ；l a r g ee d d ys i m u l a t i o n ；t u r b i n e ；b l a d e h y d r o f o i l ；h y d r o d y n a m i c sa n a l y s h “： q ： f ： f 凹： y ： f = 1 i ，： ： 瓯： ： 4 ： c ：互： 4 c c ，2 万： 屹o p = 厶： “= 吉莓五： p - ；4 ： r e ： m a 丘： 主要符号说明 与气体分子质量、动量、能量相关的函数 粒子速度矢量 粒子空间位置矢量 粒子加速度矢量 b o l t z m a n n 方程的碰撞项 粒子速度分布函数 粒子速度平衡态分布函数 碰撞频率 弛豫时间 权系数 单位格子横向长度 单位格子纵向长度 单位时间 粒子速度 格子声速 格子运动粘度 密度 速度 压强 雷诺数 马赫数 计算区域长度 t ： 虬、u o ： q ： q ： q ： l ： d ： 厶： u r ： ： a o a 、口： 占：互：旦： f dc d c ： 五： 仃： q ： c q ： k p ： p 毛p a v 3 ： 计算区域宽度 远端入口流速 压力系数 升力系数 阻力系数 升力 阻力 尺寸比例因子 速度比例因子 时间比例因子 翼型攻角 升阻比 翼型弦长 尖速比 密实度 水轮机转矩 水轮机转矩系数 水轮机获能系数 轮机轴输出功率 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 j - - - l 刖舌 1 1 研究背景 1 1 1 海流能利用现状与开发的关键技术 能源是人类社会存在和经济发展必不可少的物质条件之一，近几十年来， 随着煤炭、石油等化石能源的日益枯竭和环境污染加剧，使得环境友好、可以 永续利用的可再生能源的有效开发与利用提上日程。全球的海域资源辽阔，蕴 藏着丰富的海洋能资源，主要包括潮汐能、波浪能、海流能、海水温差能、海 洋盐差能等，其可再生能源理论储量可达7 6 6 亿千瓦【l 】，可见其资源十分丰富， 利用价值极高。 作为海洋能的重要组成部分，海流能【2 】是指流动的海水所具有的动能，主 要指海底水道和海峡中由于潮汐导致的有规律的海水流动而产生的能量。海流 能的能量正比于海流流速的平方和流量，具有能量密度大，载荷稳定，预测性 强等优点，而且全球蕴藏的具有开发利用潜能的海流能总量达3 亿千瓦【3 】。 海流发电技术是目前海流能开发利用的重要技术之一，其中海流能获取装 置将海水的动能转化为机械能，进而转化为动能，是海流能发电技术的关键部 分，关系着整个发电系统的性能。目前大部分使用的海流能换能装置为水轮机， 其水动力性能优劣直接影响了海流能有效的开发利用，是海流能开发利用的关 键水动力构件。水轮机主要由主轴、叶片以及相关的连接控制机构组成，根据 其结构运动的情况可分为水平轴和垂直轴两类【4 】，如图1 1 和1 - 2 所示： 图1 1 水平轴水轮机 基r 格子b o lt z m a n n 方法的海流能水轮机人雷诺数模拟研究 图1 - 2 垂直轴水轮机 水平轴水轮机旋转轴向与水流方向平行，当水轮机与水流方向没有夹角时， 发电功率较大，有夹角时，获能利用率下降，此类具代表的水轮机是英国m c t 公司的s e a f l o w 和s e a g e n 5 】：垂直轴水轮机分为竖轴水轮机和横轴水轮机两种， 竖轴水轮机的旋转轴向与水流方向和海平面均垂直，不受水流方向的影响，美 国东北大学g o r l o v 教授发明的g o r l o v 式水轮机1 6 是典型代表：横轴水轮机的旋 转轴轴向垂直于水流方向平行于海平面，获能利用率受水流方向影响，因其应 用的局限性故未被广泛使用。 水轮机叶片翼型的设计是海流能水轮机的关键技术之一，其水动力性能影 响着水轮机的获能效率，通过对水轮机叶片翼型进行优化，获得具有较大的升 力系数和较小的阻力系数的叶片翼型能提高水轮机的获能效率。 1 1 2 格子b o l t z m a n n 方法介绍与研究现状 近些年来，计算流体力学( c o m p u t a t i o n a lf l u i dd y n a m i c s ，简称c f d ) 的 数值模拟方法在工程中的应用越来越广泛，利用计算流体的数值分析方法对水 轮机及其翼型叶片进行模拟研究分析可以为水轮机的设计优化提供理论依据。 由于传统的数值分析方法存在计算速度慢、收敛效率低，可扩展性差，对计算 机资源要求较高等缺陷，因此亟需寻找更为有效的数值方法。 格子b o l t z m a n n 方法( l a t t i c eb o l t z m a n nm e t h o d ，简称l b m ) 诞生于2 0 世 纪8 0 年代末，格子b o l t z m a n n 方法不同于传统模拟方法不同，它基于分子动理 论，它在宏观上离散，在微观上连续，是一种介观模拟方法。b o l t z m a n n 方程 比n a v i e r - s t o k e s 方程更为基本的反映流体机理，在2 0 年的发展过程中，也已 成为一种极具发展前景的模拟手段。在非牛顿流体忉、粒子悬浮流【8 1 、湍流【9 ，1 0 1 、 化学反应流【1 1 1 等许多领域得到应用。 基于格子b o l t z m a n n 方法的海流能水轮机大雷诺数模拟研究 格子b o l t z m a n n 方法具有微观方法上假设条件较少，宏观方法上不关心分 子运动细节的优势【1 2 1 ，该方法相比较于传统方法来说具有许多独特的优势：在 算法上具有数值稳定性，有更高的数值精度；边界条件处理简单；压力项求解 简单，不需要积分运算，适合于求瞬态值；由于其网格特性适用于并行计算， 可以大大提高运算效率等等。随着格子b o l t z m a n n 方法理论的完善与发展，其 在工程实际中的应用也越来越广泛，因此探索将该方法应用于海流能水轮机的 水动力特性分析问题上具有重要的价值和意义。 1 2 国内外相关的发展与研究现状 1 2 1 国内外海流能发电技术现状与进展 国际上海流能发电技术的研究源于2 0 世纪7 0 年代中期，经过3 0 多年的发 展历程，越来越多的科研机构和公司参与到海流能发电技术的研究中，海流能 发电技术日趋成熟，目前国际上已有的海流能发电装置的水轮机主要有：英国 m c t 公司研制的水平轴海流能水轮机s e a f l o w 和s e a g e n ；l u n a re n e r g y 公司 设计的多叶片转子的水轮机r o t e e ht i d a lt u r b i n e ( 啪；爱尔兰o p e n h h y d r o 公 司开发的o p e n c e n 仃e 水轮机；英国d t i 设计的竖轴变倾角水轮机等等。 我国的海流能发电研究始于2 0 世纪7 0 年代末，由于近年来海洋可再生能 源的开发利用受到高度重视，海流能发电技术的研究与应用也取得了良好的发 展，东北示范大学、哈尔滨工程大学、浙江大学、中国海洋大学的研发团队也 相应取得了一定成果，海流能发电技术研究有望很快进入工程示范和规模化应 用阶段【1 3 1 。 1 2 2 格子b o l t z m a a n 方法在该领域的应用研究 目前，格子b o l t z r a a n n 方法在工程实践中的应用越来越多，在水力学问题 的模拟上也已取得一定成果。如f i l i p p o v a 等用局部网格加密计算圆柱绕流，在 r e = 1 0 0 时的s t r o u h a l 数为0 2 9 7 ；美国德拉维尔大学机械工程系( d e p a r t m e n to f m e c h a n i c a le n g i n e e r i n g ， u n i v e r s i t yo fd e l a w a r e ) 利用该方法实现在温度分布 传感器监测油井生命周期工程应用中弯曲管道中层流和湍流的模拟，如图1 3 所示；武汉大学的程永光用非均匀网格计算圆柱绕流，通过不同时刻的流线图 反映了漩涡产生到脱落的过程，展现了j 6 f 现象，与实验观察的现象非常一致， 同时他还将该方法应用于水利工程，提出了一维、二维明渠非恒定流模型和一 3 丛j i 格rb o lt z m a n n 办法的海流能水轮机人雷诺数模拟研究 维水利过渡过程模型，计算了工程实际中一维圆形溃坝波的传播和一维弯道水 锤，特别是一维钢筋混凝土蜗壳水锤的计算，得到了清晰的压力变化过程和流 速分布图形，实理卜海水坝、= 峡t 稃水轮机的模拟( 图1 4 ) 等 訇1 3 温度分布传感器监测油井生命周为 图l4 三峡工程水轮机( d 1 - - 9 5 m ，r e = 1 0 8 ) 在海流能水轮机发现技术研究中，叶片翼型的优化是一个关键技术，目前 对翼型的模拟分析也已取得了一定成果 1 4 1 7 1 ，如西北工业大学通过插值补充 l a t t i c eb o l t z m a n n 方法实现了复杂的低雷诺数下翼型绕流的模拟，t a r oi m m u r a 等将插值l b m 与b a l d w i n l o m a x 湍流模型结合，模拟了n a c a 0 0 1 2 在雷诺数 r e = 5 x 1 0 5 时的绕流，得出升阻系数等。 1 3 本文的主要研究内容 1 课题来源 本学位论文得到国家自然科学基金青年科学基金资助项目“基于l b m 方法 的潮流能水轮机转子水动力学特性研究( 5 0 9 0 6 0 7 5 ) ”和山东省优秀中青年科学 4 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 家科研奖励基金资助项目“潮流能水轮机水动力学特性的l b m 方法研究”的支 持，旨在解决l b m 用于海流能发电装置水轮机水动力学特性分析中l b m 对于 工程实际应用中大雷诺数下模拟的应用限制，有效解决水轮机转子复杂流场模 拟的一些难点问题，可以进一步对该方法的工程应用领域进行拓展，并通过对 水轮机的性能研究为水轮机转子优化设计提供理论依据。 2 研究内容 本研究拟采用l b m ，分别从网格划分，进出口边界处理，上下边界的边界 条件，碰撞过程以及与湍流模型相结合五个方面探讨实现大雷诺数下海流能发 电水轮机翼型叶片在海流中的力学特性的模拟分析。对不同的设计参数进行求 解，模拟求解大雷诺数下单个叶片翼型和多叶片水轮机的水动力学参数，研究 其水动力特性。 3 拟解决问题 在实际工程问题中，往往是湍流流场，要将l b m 应用于实际工程中，模 ￥ 拟大雷诺数的流场是非常必要的，针对目前l b m 的应用与研究，本文采用l b m 与湍流模型结合的方法来模拟大雷诺数下的流场；通过合适的网格划分与边界 条件对单个叶片翼型进行模拟分析，通过程序运算求解了翼型水动力参数并导 出数据进行分析，得到水轮机翼型叶片水动力参数在不同参数下的变化趋势； 对于多叶片水轮机的整体模拟，由于水轮机在水中是运动的，需实现水轮机的 动态模拟以及水动力特性分析。 4 研究方案与技术路线 ( 1 ) 研究方案 本论文从理论学习和理论应用两方面入手，其中理论学习包括l b m 基本 理论知识和翼型相关理论知识，特别是大雷诺数下绕流模拟的程序搭建；理论 应用主要指将l b m 应用于水轮机翼型叶片以及水轮机整体的水动力特性分析 上的研究，研究方案框图如1 5 所示。 5 基于格子1 3 0 lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 图1 - 5 研究方案 ( 2 ) 技术路线 具体实现步骤如下： a 翼型模型数据生成与翼型绕流网格划分与边界处理选择； b 标准l a t t i c eb o l t z m a n n 算法及程序实现，计算其所能模拟的流场雷诺数； c 湍流模型与格子b o l t z m a n n 算法结合及其程序搭建，实现模拟大雷诺数 流场的程序： d 单个水轮机叶片翼型绕流动力学分析与对比验证； e 不同叶片数量的垂直轴水轮机的绕流模拟及水动力性能分析； 利用l b m 方法探讨提高水轮机发电效率的措施。 其技术路线图如下图所示： 锵溪 实验验证了司秣_ i = i 优化 _ 一输 h 绐粜 图1 - 6 技术路线图 6 翼礴多数优化 获撄轻性援爨 蔫 氍一靖一 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 2 格子b o l t z m a n n 方法的理论与数值模型 2 1 从b o l t z m a n n 方程到格子b o l t z m a n n 方程 2 1 1b o l t z m a n n 方程 b o l t z m a n n 方程是气体动理论基本方程，是描述离子速度分布函数厂时空 变化的守恒方程【1 8 】。因为在任何一个宏观系统中，每个粒子的微观运动都遵循 着力学的基本规律，因此只需计算出大量粒子的个别运动，系统的宏观参数就 能确定；而b o l t z m a n n 方程的基本思想是不去确定每个粒子运动状态，而是通 过求出其在某一状态下的概率，进而通过统计方法确定系统的宏观参数0 9 。其 推导基于如下三个重要的假设：分子互相的碰撞只考虑二体碰撞；分子混 沌假设，认为粒子在碰撞之前不相关；外力不影响局部碰撞的动力学行为。 b o l t z m a n n 方程的表示如下： 等+ 善考+ 口差= ( 研一觑1 ) d z l g l = o s o d r l d 磊 ( 2 - 1 ) 其中，厂表示粒子的速度分布函数，孝是粒子速度矢量，是粒子空间位 置矢量，t 为时间的函数，a 为加速度矢量，右端项称为碰撞积分或碰撞项 ( c o l l i s i o nt e r m ) ，通常因为嘶) 这一积分的存在给b o l t z m a n n 方程的求解带来 了极大的困难。 b o l t z m a n n 方程与流体力学基本方程有密切联系，其中m a x w e l l 输运方程 是联系b o l t z m a n n 方程与流体力学宏观基本方程的桥梁。如果将b o l t z m a n n 方 程两侧同时乘一个与气体分子质量、动量、能量都相关的函数，然后积分即 可以推导出m a x w e l l 输运方程【1 2 1 ，而且右端碰撞项可以转化得到一个碰撞前后 物理量之差的乘积形式，此形式对于碰撞前后物理量的守恒具有重要意义，由 此可以得出：如果碰撞前后守恒( 质量、动量、能量) ，那么由这些物理量组成 的函数满足m a x w e l l 输运方程，而且右端为0 ；分别将= 朋、- - = 蟛、 = 去蟛2 代入m a x w e l l 输运方程，则可以推导出流体力学的各宏观方程，如连 续方程、动量方程、能量方程1 2 0 。 7 基于格子1 3 0 lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 2 1 2b g k 近似 b g k 近似【1 9 】是在1 9 5 4 年由b h a t n a g a r ，g r o s s 和k r o o k 提出的，为了解决 b o l t z m a n n 方程中碰撞项求解的困难，提出用一个简单的算子研来代替该碰撞 项蝴) ，算子9 应具有如下两个重要性质【2 i 】： ( 1 ) 对于碰撞不变项= m 、= 蟛、= 寺蟛2 ，g 应满足 r e , d e = 0 ( 2 2 ) ( 2 ) 由b o l t z m a n nh 定理有 i ( 1 + i n f ) j ( f 侨1 ) d e 0 ( 2 3 ) 于是求解b o l t z m a n n 方程的难点转化为寻找合适的算子骗认为碰撞的效 应在于改变使其趋于m a x w e l l 平衡态分布尸时得到最简单的算子纠1 9 1 ，引入 系数uv 是一个与分子速度无关的常数，其大小为改变率的大小和尸与厂的差 值的比例系数，同时引入算子徘 g = v 广( ，d f ( r ，孝，f ) 】( 2 - 4 ) 这样可以将b o l t z m a n n 方程简化为： 誓+ 孝芒+ 口篓= v ( 厂叼一力( 2 - 5 ) a | ja ra 芒 、 一 。 该方程称为b o l t z m a n n - b g k 方程。其中v 是单位时间内，处一个速度为喜 的分子与别的分子的碰撞总数，也叫做碰撞频率。另一个重要的参数为碰撞时 间，定义为碰撞频率的倒数： = 二 ( 2 - 6 ) 亦称为松弛时间( r e l a x a t i o nt i m e ) 或弛豫时间，是指两次碰撞的平均时 间间隔。 b g k 近似使b o l t z m a n n 方程线性化，求解更为方便，其意义为把气体的平 衡态的过程看做为一个简单的弛豫过程，主要用于处理偏离平衡状态不远的小 扰动问题，因为它的简单性所以应用广泛。 2 1 3 格子b o l t z m a n n 方程 格子b o l t z m a n n 方程可以看做是求解连续b o l t z m a n n 方程的一种离散格式， 8 基于格子1 3 0 lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 是b o l t z m a n n b g k 方程的一种特殊离散形式【2 2 1 ，包括速度离散和时空离散。粒 子的速度是连续的，不停的做着碰撞与迁移的运动，而且其速度在相空间上是 无穷维的，由于粒子的运动细节无法显著决定流体的宏观运动，因此将粒子的 速度简化为有维的速度空间，从而求解粒子的分布函数成为关键。在不考虑有 外力作用的情况下，我们将变化后粒子的分布函数看做变化前粒子的分布函数 与碰撞过程中粒子的分布函数之和，表示如下： 五( x + e o s t ，f + 西) = 厶( 而f ) + q 口( 厂( x ，f ) ) 口= o l ，g 一1 ( 2 - 7 ) 其中的碰撞过程项可以用分布函数五表示，e h j = 述方程( 2 - - 4 ) 和( 2 - 6 ) 可以进一步推出得到： q 口( 力= 一三( 厶一五q ) ( 2 - 8 ) 由此可得到不含外力项的格子b o l t z m a n n b g k 方程( l b g k ) ： 五o + 乞研，f + 衍) 一五( 五f ) = 一三眈( 五f ) 一z 叼( 茗，t ) l ( 2 9 ) “ l b g k 方程具有l a g r a n g e 特性，其数值计算具有二阶精度，由粒子的离散 速度我们可将粒子的运动在物理空间分为碰撞和迁移两个相对过程。在该方程 中，平衡态分布函数的求解是关键，不同的网格剖分会得到不同的平衡态分布 函数，所以确定不同的网格所对应的平衡态分布函数是应用l b m 建立模型的 关键。 翟 2 2 格子b o l t z m a n n 方法的基本模型 2 2 1 平衡态分布函数 一个完整的格子b o l t z m a n n 模型由三个部分组成伫3 1 ：离散速度模型即我们 通常所说的格子；平衡态分布函数；分布函数的演化方程。构造l b m 模型的 关键是选择合适的平衡态分布函数，1 9 9 2 年，q i a n y u e - h o n g 等人提出的具有d 维空间q 个离散速度的d d q q 系列模型是l b m 的基本模型。该系列采用如下 的平衡态分布函数： f f f = p 1 + 虿e a u + 簪一身 ( 2 - 1 0 ) 其中，最常用的二维模型和三维模型分布为d 2 q 9 模型和d 3 q 1 5 模型。 d 2 q 9 即二维九点正方形网格模型，如图2 一l 所示，该模型的离散速度如 9 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 下： 气2 ( o 嘭口= 0 唯醐p 峥，啡一1 ) 2 1 , 2 , 3 , 4 ( 2 - 11 ) 风绯一呼，s i n 【f 2 a 一峥弦一5 , 6 , 7 ，g 可得其对应的平衡态方程为： h 4 + 学+ 掣一纠 协2 ， 其中，w o = 4 9 ，m = k = 9 ，w 5 = w 3 = 1 3 6 。 d 3 q 1 5 为三维十五点的网格模型，如图2 2 所示，对用的离散速度向量矩 l0l 一10 0o ol 一1111一l11l e = 10 0 0 1 1o o1 一l 11 11一ll l ( 2 1 3 ) 10 00001 11 1 一l l11 一l1l 其中，权系数和声速分别为：= 吾，q = 万1 ，鸭= 瓦1 ，q = 万1 。 i 7 x 图2 - 2 d 3 q 1 5 网格模型 通过c h a p m a n - e n s k o g 展开可以导出每个格点上的宏观变量对应的方程， l o 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 模型的宏观密度与速度可以通过粒子分布函数得到，其定义如下： 户= 军石，“= 吉军石q 协1 4 ) 格点上的压强通过密度来求解： p = 譬c 2 = 硝 ( 2 1 5 ) j 其中，c j2 去为声速。 弛豫时间与运动粘性有关 l ，：2 r - 1c 2 a t ( 2 1 0 一) l ，= 一c 。 l z lj 6 其中，1 ，为运动粘性。因为已经选定了合适的格子尺寸和离散速度，那么 运动粘性1 ，可以根据给定的雷诺数来计算。雷诺数、运动粘性之间的关系如下： r e ：u o l( 2 1 7 ) 其中，砜为远场速度，l 为特征长度。弛豫时间r 可以根据公式( 2 1 6 ) 求出： f ：二坚+ 0 5( 2 1 8 ) r e c 。a t 2 2 2 网格技术 格子b o l t z r a a n n 方法建立模型时对流场的网格划分有两种方式：均匀网格 和非均匀网格。 均匀网格是标准格子b o l t z m a n n 方法主要采用的网格划分方式，通过碰撞 迁移步骤进行求解，是格子b o l t z m a n n 方法领域应用最为普遍的一种网格技术。 均匀网格划分即离散速度等于格子间距与时间步长的比值，从而保证粒子碰撞 迁移之后会落在下一个网格节点上【1 2 1 ，该网格划分方法简单但要求网格对阵划 分而且为了获得局部区域精确解，需要加密整块计算区域网格，但是这样会导 致所需硬件资源与计算时间的提高，限制了l b m 方法的应用。 由于均匀网格划分技术的缺陷，非均匀网格技术应运而生，获得非均匀网 基r 格子b o lt z m a n n 方法的海流能水轮机人雷诺数模拟研究 格的方法有两种【1 4 】： 一是在需要高精度求解的局部区域上建立较好的细化的网格。在l b m 中 采用自适应技术( a m r ) 2 4 1 或者复合网格系统【2 5 1 。虽然网格之间的间隔不均匀， 但每个网格的运算法则与l b g k 是完全相同的，只需在边界网格处进行特殊的 处理。 二是采用贴体网格( 见图2 3 ) ，可以通过插值方法、泰勒展开、有限差分 等方法获得。详细推导过程可以参照何雅玲教授等人编写的书籍【1 9 】。 图2 - 3l b m 中的非均匀网格：贴体网格( 】4 j 2 2 3 边界处理 边界处理在求解物理问题中非常重要，对于同样的问题，由于计算方法和 模型不同，所使用的边界处理方法也不同【1 4 】。采用l b m 进行计算时，流场内 部的分布函数可以通过分布函数的演化方程求得，但边界节点上的部分分布函 数是未知的，例如在出入口边界，流固耦合的边界位置等，需要根据已知的宏 观边界条件确定边界节点上的分布函数才能进行下一步的计算，求解这一过程 的方法即为格子b o l t z m a n n 方法的边界处理方法。l b m 适用于很多边界条件的 处理方法，可以容易的处理复杂的边界条件。 在l b m 方法中，边界条件可以分为速度边界和压力边别1 9 】，边界条件划 分如图2 4 所示，下面介绍几种常用的边界处理的方法。 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 图2 4 边界处理格式 1 反弹边界 反弹边界( b o u n c e - b a c k ) 是启发式格式中的一种，也是格子方法中处理固 壁边界时最常用的办法，包括标准反弹格式、半步长反弹格式、修正反弹格式 三种。该边界处理的基本思想是保证壁面处的流体宏观速度为o ，即对静止的 固体边界上的粒子作弹回处理。 粒子在碰撞、迁移的过程中，每一次碰撞后从流体迁移向壁面的粒子是已 知的，而从壁面迁移向流体的粒子是未知的。假想从壁面迁移出的粒子就是向 壁面运动的粒子的反弹，图2 5 描述了d 2 q 9 模型中的反弹边界处理示意图。 n 。l 哌 夕 图2 - 5d 2 q 9 模型中的反弹边界处理 由于反弹边界在计算中非常容易实现使得l b m 在处理复杂几何边界问题 时具有非常大的优势。 2 非平衡态外推格式 非平衡态外推格式是外推格式的一种，由g u o 等人 2 6 1 提出的。该格式的指 导思想是将边界格点上未知分布函数分解为平衡态和非平衡态两部分。平衡态 部分由边界条件的定义获得，同时利用一阶精度的外推格式得到非平衡态部分。 该格式可以表示为 刀一刀叼= 石1 一彳凹 ( 2 1 9 ) 1 3 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 其中，o 和z 1 分别代表边界上的格点与其紧邻流体格点处的分布函数， ，o 。钾和爿凹分别为其相对应的平衡态分布。 3 z o u h e 边界处理方法【明 ( 1 ) 通用性z o u h e 边界条件( 已知畋和”，) z o u 与h e 在非平衡态反弹边界条件基础上提出了一种新的边界条件：如 图2 - 6 所示，比如在下边界，经过迁移之后石，z ，石，五，石，五已知，假设在壁面 蚝和已知，则可以通过求解宏观密度和速度的公式求出五，石，五，其计算过 程如下： a 、 i 。 ，7 、r 。e 图2 - 6z o u - h e 下边界处理 五+ 石+ 五= p 一( 石+ 石+ 石+ 厶+ 石+ 五) 正一五= 肛一( 石一石一石+ 石) 五+ 五+ 五= 以+ ( 厶+ 石+ 五) 由以上方程求得密度表达式： p：1【石+石+石+2(五十石+五)】l p 2 - u , 【厶+ 石+ 石+ 2 u 4 + 石+ ) 】 ( 2 2 0 ) ( 2 2 1 ) ( 2 2 2 ) ( 2 2 3 ) 假设反弹边界条件仍然适用于非平衡态部分的分布五一正叼= 石一五吲， 由此求得正，石，五： 五= 厶+ 詈以 石= 石一三研一石) + 互1 纵+ i 1 鹏 1 4 ( 2 2 4 ) ( 2 2 5 ) 基于格子b o lt z m a n n 方法的海流能水轮机大雷诺数模拟研究 兀= 石+ 圭( 石一石) 一i i + 否l 鹏 ( 2 2 6 ) ( 2 ) z o u h e 压力边界条件 z o u h e 压力边界条件采用与速度边界条件类似的方式处理压力边界条件。 以入口边界为例，如图2 - 7 所示，设流体边界沿着垂直轴】，方向，压力已知( 密 度已知p = 彬) 。b 已知为0 ，= o ，迁移处理后，五，石，五，以，石已知，需 要确定石，石，石 z + 石+ 石= 几一c o + 厶+ 石+ 六+ 五+ 石) ( 2 2 7 ) 石+ 石+ 石= 风略+ c a + 五十石) ( 2 2 8 ) 石一五= 五+ 工一五+ 石 ( 2 2 9 ) 然后再利用z 一彳q = 石一z 凹得到如下关系式： 略：l 一垃堕丛兰垃趔( 2 3 0 ) 心 从而可以求得： 石= 石+ 成 ( 2 - 3 1 ) 石= 石一圭( 正一工) + i 1 忍致( 2 - 3 2 ) 石= 五+ 丢( 五一石) + 吾1 以叱 ( 2 3 3 ) _ z 红_ z 。芒。己z z 丝z z z 2 己己z 己z z z i 己z z 岔 。0 酶l l 旺 啦 图2 7 入口边界 同理对于出1 2 1 边界，如图2 - 8 所示，经过迁移后 ，五，工，石，石已知，需 要确定；，五，石，应用上述原理及推导方法可以求得： 1 5 基于格子b o l t z m a n n 方法的海流能水轮机大雷诺数模拟研究 z = 彳一岛心 ( 2 - 3 4 ) 以= 石一圭( 五一工) 一i 1 岛叱 ( 2 3 5 ) 石= 五+ 三( 五一工) 一否i 岛( 2 - 3 6 ) 在本论文中所进行的绕流模拟的边界条件均是采用z o u - h e 压力边界条件。 - - - l - _ _ 呻 3 秘l 耗t 图2 - 8 出口边界 2 3 格子b o l t z m a n n 方法与实际流场之间的单位转换 将l b m 应用与工程实践的一个重要前提是解决格子单位与物理单位之间 的转换关系，在应用l b m 方法进行数值建模与分析时一般采用的都是格子单 位，为使数值模拟时的模拟参数更加接近于实际工程应用，了解格子单位与实 际流场之间的单位转换是关键，通过两者之间的转换关系，可以在进行数值模 拟的时候选择与工程中相对应的模拟尺寸进行研究分析，s u c e i 教授2 8 1 与何雅 玲【1 9 】教授在这个方面都已进行了一些探讨与研究。 格子单位与物理单位之间的转换的基本思想是确保无量纲化前后的流动与 换热准则数一致，也即根据流体力学的相似性原理，在保证流场的关键准则数 如r e y n o l d s 数相同的情况下，可以将流场的性质看做是相同的，在l b m 模拟 计算时也是这样的基本思想，保证模拟的关键数准则与实际工程问题中的流场 的相关准则数相同，在进行水轮机的模拟中，需要保证r e y n o l d s 数的不变。此 外，在格子模拟中的另一关键参数弛豫时间【2 9 捌r 的表示为： f ：皂+ 0 5 u( 2 一- 3 一7 )f = + lj 其中，格子声速c 。= ；，粒子速度c = 睾，c 确定之后可以确定格子中粒 。4 3 访 1 6 兰一 基于格子b o l t z m a n n 方法的海流能水轮机大雷诺数模拟研究 子每个方向的运动速度。 从而可以得到格子运动粘度屹 ： 屹= 竽c = 竽鲁 协3 8 ， 屹2 i c 2 了一吉 ( 2 - ) ，1 在实际工程应用中，可知海中的实际物理粘度与声速分别为： y = 1 0 x 1 0 。6 m 2 s ，e = 1 5 0 0 m s ， 由蚱= 曼可以求出速度比例因子【1 2 1 。 乞 由蚱= 二，可以得到尺寸比例因子 屹 争铥v - - 2 毒= 旦c c ：( r - o 5 ) 像3 ” 、1