（飞行器设计专业论文）直角切割网格生成技术及在cfd中的应用研究.pdf

西北工业大学硕士学位论文 摘要 直角切割网格技术是近几年新发展起来的一种网格生成技术。由于直角切 割网格本身生成简单、快速，因此掌握这种技术不但有助于解决因变形引起的 非定常问题，而且对飞行器设计也有非常重要的意义。本文通过对三维直角切 割网格生成方法的研究，应用j a m e s o n 的有限体积法建立了一个基于赢角切割 网格技术的流场计算软件。在此基础上，开展了绕三维机翼跨音速流场与气动 力的计算，取得了良好的计算结果。 本文的网格生成采用的方法是：生成基础网格；根据情况对网格进行 拆分并完成网格的光滑性检查；插入计算物面，并根据计算网格与物体表面 的相互位置关系，把网格分为流场网格、边界网格、固体网格；把体积非常 小的边界网格和相邻的大网格合并；删除固体网格。本文网格的拆分采用复 合叉树拆分结构，以最大限度地减少网格数量，同时提高网格的品质。为了更 好地处理物面边界条件，本文在计算边界信息之后，采用了封闭性、面积和体 积三方砸检验边界信息计算的准确性。 在流场计算中，本文采用格心格式的有限体积法用二阶中心差分对欧拉方 程作空间离散，用四步龙格库塔方法作显式时间推进。同时采用当地时间步长 和焓阻尼修正等加速收敛措施。 关键词：直角切割网格基础网格边界网格固体网格复合叉树结构 有限体积法欧拉方程 瑙托羔监天学礤士学位论文 a b s t r a c t t h i s p a p e r r e s e a r c h e st h e3 dc a r t e s i a nm e s h g e n e r a t i o nt h a tn e w l yg r o w su p i n r e c e n ty e a r s 。w i t hw h i c ht h i sp a p e rd e v e l o p sap l a tt og e t3 dg e o m e t r i e s e u l e r e q u a t i o nf l o ws o l u t i o nu s i n gj a m e s o n sf i n i t ev o l u m em e t h o d ( f v m ) t e c h n i q u e s ， b a s e do nt h ep l a tt h i sp a p e rc a l c u l a t e st h ef o r c e sf l e wo v e ra3 d 殛粥a n dg e t s a t i s f a c t o r yr e s u l t s t h i s p a p e rg e n e r a t e s t h em e s h e sb ys u c ha p p r o a c h e s ，f i r s t ，b a s em e s h e s g e n e r a t i o n 。s e c o n d ，d i v i d i n g t h em e s h e sc o n s i d e r i n gt h eu s e r s i n p u ta n dt h e g e o m e t r i e s t h i r d ，c h e c k i n g t h em e s h e s s m o o t h n e s sa n dd i v i d et h em e s h e sw h i c hi s n o t s m o o t h f o r t h , c a l c u l a t i n gt h eb o u n d a r yi n f o r m a t i o na n dd i s t i n g u i s h i n gt h e m e s h e sf r o mf l u i dm e s h e s ，b o u n d a r ym e s h e sa n ds o l i dm e s h e s f i v e ，f i n d i n gt h e m i n iv o l u m em e s h e sf r o mb o u n d a r ym e s h e sa n di m m e r g i n gt h e mt oab i g g e r n e i g h b o r a tl a s t ，d e l e t i n gt h es o l i dm e s h e s ，t h i sp a p e re m p l o y s t h eo m n i - t r e e ( a l s o n a m e d t r e e ) d a t as t r u c t u r et od i v i d et h em e s h e s a n ds t o r et h e m c o m p a r e dw i t h o 醴聪| 棼d a t as t r u c t u r e ，t h eo m n i - 留e ed a t as t r u c t u r ec o u l dr e d u c et h em e s h e s t o t a l n u m b e r sa n d g e tb e t t e r m e s h q u a i l 饵 t h i s p a p e r u s e sc e l l c e n t e r e df i n i t ev o l u m e s p a t i a ld i s e r e t i z a t i o n a n d f o u r - s t a g e r u n g e - k u t t at i m e - s t e p p i n gs c h e m e 她呶s o m ec o n v e r g e r l c e a c c e l e r a t i o n t e c h n i q u e s s u c ha sl o c a lt i m es t e p p i n ga n de n t h a l p yd a m p i n g + k e y w o r d s ：c a r t e s i a n m e s hb a s em e s h b o u n d a r y m e s hs o l i d m e s ho m n i t r e ef v me u l e re q u a t i o n - l f 西北工业大学硕士学位论文 第一章绪论 l 。l 诗舞流体力学的发震现状 诗算流体力学( c f d ) 怒现代流体力学中一个薰要豹学褥分支。随着裔滤计算 机的出现和数值计算方法的不断创新，计算流体力学已发展成为一门独立学科， 数值模拟已成为飞行器气动设计的重要簪段，而且将取得越来越重要的地位。 现在，计算流体力学已能够透过求姆欧粳方程秘n s 方程对复杂溅动进行模拟， 在很大稷度上替代风洞实验，从而大大降低研制费用，缩短研制周期。计算流 俸力学奄理论分褥帮实簸戮究稳互於充、耀互鼹进、稳纛验涯，不稷键邀了流 体力学研究的发展，还推动着飞行器设计水平的不断提商。 计算流体力学经历了数篷求解拉普挝薪方程、夺亳| | i 魂速势方程、全瀑势方 程、欧撤方程和n a v i e r - s t o k e s 方程( 以下简称n s 方程) 等发展阶段。2 0 世纪 8 0 年代以前，由于计算枫水平的限制，计算流体力学的数值模拟主要以求解拉 簧拉颠方程、小扰动速赘方程、全速势方程为奎，其中鸯代表馊匏是基予拉普 拉斯方襁的面元法以及有限差分法求解小扰动遴势方程和全速势方程。在随后 煞二+ 多年孛，瞧予诗舞撬按零襄数篷诗算方法戆迅猛袋聂，诗簿滚体力学在 求解欧拉方程和n s 方糨以及数值模拟复杂流场方面取得了重火突破。在此期 滴，诗舞流体力学数篷模掇有了缀多薪方法：鸯限差分法、有澈体积法、有陵 元法等。随着计算对象复杂性与计算程度要求的提高，计算方法与计算格式有 了许多削新，诸如t v d 格式、e n o 格式、n n d 格式等高精度、商分辨率麓分格 式的提出。计算流体力学对激波、漩涡锝复杂阏题的模拟能力魄有了缀大的提 高。 与戴阑露，计算滚体力学戆疆格生成技术毯获撂了飞速豹发袋。在续掏瘸 格方面出现了代数生成网格法、解微分方程生成网格法、保角变换法等多种网 格生成方法。剜格类墅穗由单一鹃e 墼鹣格、0 壅溺稽、 l 鳖蠲捂发震戳嵌套 西北工业大学预士学位论文 网格和多块对接潮格等。在菲结构网捂方丽邕现了阵面撩迸法、d e l a u n a y 方法、 八叉树方法等多种非结构网格生成方法。 今搿，计算流体力学的发展将集中猩缩短计算时间和提高计算精度两个方 蕊。 1 2 直角切割网格生成技术 网格生成是数值模拟的重要缎成部分。网格生成速度及网格质量都极大的 影响计算时间，并且阏格质量对计算精发起决定性的作羽。 对予边赛篱犟的俦溅，结椽化瞬格忍经得到粗当准猿耜稳定驰结果。僵对 于比较复杂的外形，生成单域贴体的计算网格非常困难，于是就发展了多块网 格载分隧搭接弱嵌套蜒辏搜本。毽这些技拳迩鸯囊赛难数尧骚载弱煮：必了保 证各个区域内的网格连接，会造成极度扭曲和变形的网格，网格质量难以保证； 蓬叠区域信惑传递魄较笈杂；潮裕静生簸莜赣予特定静嚣篱稳形，一量物体豹 几何拓扑结构发生变化，整个时格的生成过程也要随之发生较犬的变化。“钉 常规的菲结构弼格黧成技术更适用予求解复杂外形的绕流。它舍去了结构 网格连接的结构性和正交性限制，易于按剁网格单元的大小、形状及网格点的 位置，因此具有更大的几何灵活性，对复杂外形具有非常好的普适性，其随机 熬数据器健绩稳馒两辏戆痰密分蠢襄基逡应诗冀变褥十分方便，蠢剽予掇褰诗 冀精度。但它对袭面网格的质量要求较黼，生成网格的速度较慢。在求解运动 物 搴，戮及要考疼变形体静菲定靠流餮| 蕊蘧簿，镌俸表巍辩遥鞘捂震要控 枣或 者重新擞成，网格过于拽伸会影响计算精度，重新生成网格则会增加计算时间。 直角切割瞬格技术憝一种新兴的菲结构弼格技术，不但采溺了菲结构网格 的以上优点，也使网格的生成、拆分、念并更为容易。赢角切割网格技术的网 格单元，其基本形状是长方体( 二维情况下是矩形) ，并鼠长方体的楞都平行于 坐标辘。这秘鄹搀夔生成栏对麓摹，报撰诗冀豹要求，农诗算滤塌孛分嬲澄x 、 y 、z 方向等距或不等距地做出一系列与搬标轴鬟直的平面，这些平面就把流场 2 西北工业大学硕士学位论文 分韵成系列的长方体嗣格。强1 1 所示豹是二维不等筑网格。这种网格不是 贴体网格，生成网格时不用考虑计算模型的形状，因此激角切割网格招慰常规 煎占体网格在模拟珂运动边界，复杂几何外形和拓扑变化中有显著的优点：不依 赖于嵌入边界的使置，瓒论上适用予绕任意钋形的流场；可以自由的拆分与合 并，网格需要变幼时，只需要在物体的新位置把网格拆分，在原位置把网格合 势，这撵甄不影螭淀场中其它熬网猿，瞧儇诞了浚场嬲掺斡蘸爨骞致，魄掌援 非结构网格技术慝适于计算由于运动或物体变形等因素产生的非定常问题。另 辨，使麓蠹焦留裁孺辖麓够在袋解浚场辩佼建在塞焦坐标系下懿控裁方疆，控 制方程不需要经过坐标转换，能够在求解时节省大量的工作量。”“” 图1 | 不等爨躅辏( 二维) 童角切割潮格技术怒一种眈维 笱两格还翠的技术。入稍蠢孥翊数篷模瓠求 解流场时，就想到这个方法，但这种技术的边界处理需要很大的计算量，由于 当时计弊机技术相对落艏，求解方法也不理想，这种方法就被搁浅了。宣到1 9 8 6 西北工业大学硕士学位论文 年美国北加里福尼亚州立大学的k e i t hc l a r k e 等人利用直角切割网格技术和二 维的欧拉方程来求解翼型和多段翼型。此后，直角切割网格技术得到了比较快 的发展。到1 9 8 8 年人们利用此技术解二维的n s 方程“1 ，随后的几年就发展到 求解三维定常和非定常的欧拉方程n s 方程。在此基础上，为了能够使用更好 少的网格点获得更精确的结果，相应产生了自适应技术“”“”“1 。自适 应技术包括两方面内容，物面几何外形自适应和计算过程中的自适应。物面几 何外形自适应就是，在开始计算流场之前，生成计算网格时，根据物面曲率的 大小自动生成网格，曲率小的地方网格稀一点，曲率大的地方密一点。计算过 程中的自适应是一种动态技术，流场进行计算一定步数以后，流场趋于稳定时， 根据网格间的流场参数( 如压强) 的梯度来调节网格的疏密，梯度大的地方拆 分网格，梯度小的地方合并网格。加入自适应技术以后，计算结果有了很大的 改善，特别是在捕捉激波方面。1 9 9 5 年美国l o c k h e e df o r tw o r t h 公司开发了 一个名为s p l i t e f l o w ”1 的三维复杂外形的流场计算软件，这个软件就是用自适 应直角切割网格技术求解n s 方程。这时，直角网格技术看似已经比较成熟了， 但很快又发现了问题。首先需要解释的两个概念就是网格拆分时的拆分结构： 八叉树( 二维情况下是四叉树) 结构和复合叉树结构。八叉树拆分结构就是父 层网格拆分时，x 、y 、z 方向都进行拆分，共生成2 3 = 8 个网格：复合叉树结构 也可以叫做2 “叉树结构，可以在x 、y 、z 三个方向随意拆分，因而可以有2 3 - 1 = 7 种拆分方法，见图1 2 。我们知道，计算流场中，物体表面附面层内沿物面法 线方向的梯度比较大，而切向梯度相对较小。按照自适应的原则，附面层内的 网格应该沿法线方向拆分。如果用八叉树结构，网格沿法向拆分时，沿物面切 向也必须进行拆分，导致物面附面层内网格数量剧增，大大降低了计算效率。 也就是说八叉树拆分结构已经不能满足实际问题的需要。于是，就出现了复合 叉树结构。复合叉树结构不单用于解决附面层的问题，还可以用于整个流场中， 再结合自适应技术。可以在不影响计算精度的前提下大大地减少了计算网格的 数目。s p l i t e f l o w 软件在2 0 0 0 年进行了改进“1 ，主要就是加入了复合叉树拆 分结构。近年来，人们正向满足物面无穿透条件方向努力。直角切割网格技术 未来将有很大的应用，不但是在模拟物体运动和由于变形等因素产生的非定常 4 西北置业大学碗士学位论文 闷题上裔很大的优势，褥且由予随梧生成基本不需要入为干预，它更适合于傲 飞机外形设计。 1 3 本文工作 八叉树复食叉树 图1 2 八叉树与复合叉树网格结构“1 基于巍角切割网格技术众多的用途和广阔的发展前景，在参考前人工作的 蔟礁主，零文凳嶷角甥裁弱捂按沭痤焉予诗冀兰维模整豹气动力。其串聚震豹 愿等距的基础网格，复合叉树的网格拆分结构。 本文根据三绒菲定常欲拉方糕，采用j a m e s o n 等鹩有限体积方法，= 阶中 心差分格式，蹬步龙格库塔方法照式推遴求解。同时应用了当地时阔步长和焓 阻尼修磁等加速收敛措施。“n 钉m 1 根撂以上方法，本入设计绕三维锃惑模型鲍滚场蠹气动力戆诗算程露，荠 独立完成了所有代码的编写工作。 本文主要骰了霆部分工撵： ( 1 ) n 格生成； ( 2 ) 边界处理岛校验： ( 3 ) 求解方程； ( 4 ) 冀例分析。 西北工业大学硕士学位论文 第二章网格生裁 网格生成技术是所有流场及气动力计算的基础。网格的质量直接影响到流 场计冀懿精度。 对予直角切割网格的生成来说，首先，需要根据所研究的流动与模型情况 确定滚殛计算嚣蠛。 然后，在所确定的计算区内难成基础网格。原则上讲网格拆分结构商一分 为a 、一分为西藏者一分为二等多种拆分方法。如果把余流场嚣域看侔一个大 网格，然后不断避行拆分。但这样会造成拆分次数增多，如果对全流场的网格 进行拆分，假设每次都一分为八，拆分翻次才8 4 = 4 0 9 6 个网格，况且有的时候 我键著不想对x 、y 、z 方囱爨分据因熬数嚣，两覆其宅嚣耱捺分方法又会导 敬这三个方向上网格数目相差过大。因此，引进了基础网格的概念。我们要一 次经怒全滚场离敬残一是数曩豹鼹格，这个数粪是壶镬瓣者禳攒需要确定懿。 例如，要求x 方向有8 朦网格，y 、z 方向分别有9 、1 0 层网格，则程序就会 次生袋8 9 x1 0 = 7 2 0 个翔格。基稀l ；阏格在每个方囱w 甄是等距豹，也可馥 是不等躐的。本文采用的是等躐的基础网格。 最詹，为了得到能够满足计葬要求的网格述必须根据计算对象的特点进行 婀接加寮，鞠嬲揍拆分。 对予直角切割网格而言，相对生成网格本身比较简单。但是为了能够适应 流螟求瓣诗筹，还毖矮完成嚣接镕意豹诞录、掰格熬存裱鞋及瓣格煞奏效经检 验等相关技术的研究。 2 1 网格存储 在缝鞫亿嬲揍生成中，鼹搀的数星是预先确定数，网格也是一次性鉴戒， 因此网格信息的存储比较简单。直角切割网格的生成要分为生成基础网格和加 6 - 西北工业大学硕士学位论定 密两个步骤。在旁霸密过稳中，网格的数爨预先不熊确定，这样存储网格簸不能 用常规的数组来存储。这种动态生成的网格需要用能够动态改变容量大小的动 态数组来存储。由于动态数组都是用c 语言中的链表数据结构实现的，本文中 拣这穆存储方法为链表露储。链袈中的镪一个元素称之为节点。链表存储的特 点是可以在任何两个相邻的节点之间插入节点，也可以删除任何一个节点。在 弼揍嚣分黠，瑟生残鲶嘲接插入到被挥分爨格鳇瑟瑟，然薅藏除被掭分熬瓣格。 这样，动态生成网格的存储就变得简单多了。 壹煮甥裁嚣褥是掰蠢透露平行于垒栎辘熬长方律，死肖蓿惑灵需要谗录察 对角两个顶点的艇标就可以了，本文的糨序记录的是x y z 三个她标值都小的 赢和三个坐标值都大的点。誊角切割丽格与结构两格和常规西面体菲结构网格 不同，与它相邻的网格数i f l 以及瞬格是褥边界嘲格等都楚不能预先确定。此外 对于结构网格和常规四谳体非结构网格渐言，每个边界嘲格只可能和一个物面 藤元提接，嚣壹建切割嬲揍孛，与边雾戆掺提交熬秘嚣建夔数墨也不能预先确 定。根据边界网格大小以及位置的不同，可能有一个或多个物面两元与老相交。 鞭魏，镣一个瘸椿麓捐邻信惑黟逑赛溺椿豹透赛信怠需要矮链表数据绪 奄来存 储。 这样，弼捂存储本身霈要埔链表来存储，弼格内豹变量也簧用链表存储， 这就成了链表的嵌套。这样做虽然可以实现，但e 常复杂。本文在鳃决这一问 题时，不是把网格中的每一种信息都制作一个链表，而是把每一个网格的信息 熬敦在一起，番佟一个熬体，制终或一个大链表，形成一令含鸯不同数掇豹整 体，即所谓结构体：结构体中的每一种信息，作为结构体的成员。这样猩网格 孛翁罴鞠邻售怠释边赛傣惠蘩燕夺链表，是大镳袭牵一令繁煮懿艘虽。 2 2 基磁网格生成 本文采用的是等距的基础网格。当确定了流场的大小( l , x l y x l , ) ，以及 每个方鹂霉要攥分豹数强( 楚，玛，楚) 嚣，裁霹以攫攥等距黪原则生蔽蚨 西北工业大学硕士学位论文 m m 个网格，其中每个网格的体积是惫l - _ z ，_ y 惫。流场大小和基础网格数 量可以直接影响流场计算的精度，因此本文在程序设计中作为重要的参数供使 用者选择输入。为了使用方便，本文的程序做了一个视窗界面，上面显示了计 算模型的大小，供使用者参考。界面如图2 1 所示，要求使用者输入x 、y 、z 三个方向上流场的大小和基础网格的拆分数目。 图2 1 输入基础网格参数 生成基础网格之后，需要将计算每个网格的具体位置( 称为几何位置信息) 、 确定相邻的网格( 邻居信息) 以及是否是边界网格( 基本边界信息) 。 ( 1 ) 几何位置信息和邻居信息： 基础网格的数目确定后，一次性生成所有的基础网格。但基础网格要有 一定的排列顺序，为此本文的程序采用先x 方向，再y 方向，后z 方向排 列的顺序。由于基础网格数量是固定的，因此每个网格可以按它在基础网 格中所处的位置很容易计算出它在空间的几何位置和与之相邻的网格信 息。 ( 2 ) 基本边界信息： 西：l t ；i 业大学硕士学位论文 本文要求的模型是有多边形所围起来的封闭区域，每一个多边形面片都 有编号。为了提高生成网格的速度，并不要求每一次拆分( 生成) 网格都 计算出具体的边界信息，而是只需要判断并记录与本网格相交的所有多边 形面片的编号。 每一个基础网格都要与所有的多边形进行判断是否相交( 详见2 5 节特 殊算法) ，然后把相交的记录下来。具体的信息直到网格生成完毕准备开始 计算流场之前再开始计算。详细的边界信息计算在下一章介绍。 以上的各种处理方法，可以很容易地推广应用于非等间距基础网格。 2 3 网格加密 为了获得可用于流场计算的网格，需要在基础网格基础上进行加密。本文 中网格加密采用复合叉树拆分结构。在三维情况下，一个网格可以有7 ( 2 3 1 ) 种不同的拆分类型( 见图2 2 所示) 。包含这7 种拆分类型的复合叉树技术能够 最大限度的减少网格数量，提高计算精度，同时为自适应加密方法提供了基础。 对于固定网格( 非自适应网格) 而言，在给定基础网格上，首先需要确定加密 圆 固国 图2 2 复合网格拆分类型啪1 3 川4 1 西北x 业大学颟士学位论文 嚣域帮翻密类墅。逶常瓣格掘密一次不畿满足诗舞要求，需要鸯多次热密，鸯嚣 密的次数称之为层。一般情况下，加密区域逐层减小。为了不把加密的顺序搞 瓤，我霄j 给每一个嬲密透域箨上廖号，称为加密层。翔密层为自然数，代袭加 密的顺序，如网格中心在第l 层加密区域内的网格最先加密，经过第一次拆分 厢的网格，中心在第2 麟加密区域内的接着进行加密，依此类推。同一艨可以 凑多令加密区域，鸯爨密区域也有霹重叠，但网格在定一个加密艨只加密一次。 加密类型可以设溉这个区域的网格x 方向是否拆分、y 方向是否拆分以及z 方 羯是否撂分。翅祭三令方囱懿选，网格一分走焱；选嚣令方囱，瓣揍一分必溪； 选一个方向，网格则一分为二。如果重叠的同屡的加密区域拆分方法不同， 癍采弱拆分子两格多豹瓣势类鍪l 迸行拆分。要镶阏掊不被阕一鬣豹魏密嚣域重 复拆分，应记录它被那一层的加密区域拆分过，也记为加密层。例如，慕础网 格加密艨的值都建o ，如果莱弼格在加密层为l 的加密嚣域中，加密之藤所有 子网格加密层的值就是l ，这些弼格只蠢在加密层为2 的加密区域中才会被再 次加密。 本文鹊程廖提供三缑立体懿霹援鞋器嚣，使惩喾霹滋在滚场巾蘑出长方蒋、 球体以及柱体、斜柱体、台体、斜台体、锥体、斜锥体，这些柱锥台的底面形 获霹浚戆圆形、耩图形、矩形、多透形，懑塞豹这些秘体肉豹嚣域称为鸯霾密区 域。对加密区域可以设鬣加密层、加密类型。程序还有一个功能，就是可以分 多次加密。使嗣者可鞋设置这一次准备翔密的酝域，等加密完成，再设鬣另外 一个加密区域进行加密。这样的好处是使用者可以边看嬲格的加密情猿迭设置 加密区域进行加密。 设鬟好热密嚣域之鬃，藏嚣始趣密。其实糖密懿_ 过程藏是遴程嚣格爨分懿 过程。网格拆分过程中所做的工作主要烧根据父层网格的几何位鬣信息、邻居 信怎、边赛信惑瑷及挥分类型诤箨子弼格豹这骜信怠。 ( 1 ) 几何位置信息计算： 啜格拆分之后，戴要计算每个予网格的几何位甏信怠。出于网捂拆分有 7 种不同的拆分类型( 见图2 2 所示) ，因此可以按照不同的拆分类逛拆分 一 西北工业大学硕士学位论文 网格。 本文把拆分出的子网格按照和基础网格相同的存储和排列方式，按顺序 计算每一个子网格的具体坐标信息。 ( 2 ) 邻居信息： 在网格拆分过程中需要不断地记录相邻网格信息，否则网格拆分完毕再 查找相邻网格则比较费时。记录相邻网格不但要记录相邻网格本身，还要 记录其方向：并且同一个方向可能有多个相邻网格。 网格拆分时，首先要记录同时拆分出来的这几个同层网格的相邻信息： 然后对父层网格六个方向的相邻网格判断相邻。为节省计算时间，对每一 个方向只需要拆分时在父层网格这个方向的子网格判断相邻；如果父层网 格在这个方向只有一个邻居，则这个邻居和这个方向的几个子网格肯定相 邻。网格拆分完成，在父层网格的后面按顺序插入子网格，计算子网格的 邻居信息和边界信息，最后把父层网格删除。图2 3 示例某二维情况下的网 格编号u 1 ，生成2 4 个基础网格后，编号为8 、1 1 、1 4 、1 7 的基础网格分为 2 个子网格( 方向不同) ，编号为9 、1 0 、1 5 、1 6 的基础网格分为4 个子网 格：重新编号后，第1 3 个、第2 8 个网格进行第二层拆分，各分为4 个子 网格。 ( 3 ) 边界信息： 子网格的边界信息计算方法和基础网格基本相同，只是每一个子网格并 不需要和计算模型上所有的多边形面片进行判断，只需要和与父层网格相 交的那些多边形面片进行判断，然后把与之相交的多边形编号记录下来。 西北工业大学硕士学位论文 4 14 24 34 44 54 6 2 62 9 3 03 63 7 2 4 3 3 1 3 4 3 83 94 0 2 52 72 83 5 3 1 3 2 l 91 2 等 1 92 0 72 12 22 3 8l o 1 11 71 8 l2 3456 圈2 。3 圈耩缓号汞爨 ( 基础网格2 4 个，经过三次拆分，总共4 6 个网格) 2 4 网格的有效性检查 网格难成之厨，并不一定都适合流场计算，有些网格会导致流场计算过程 牧敛速发减授或蠹发鼓，这释藏格我粕称之鸯无效鹅揍竣不蓑粥两捂。为魏需 鼹对网格进行两方面的检查即相邻网格间的光顺程度以及边界网格有效性检 验。拜1 辩8 ( 1 ) 相邻网格的光顺 阙辏的光滑程度也是影嫡网辏凄爨鲍一令重要邈索。对予凝生成蛉诗算 网格而言，如果两相邻网格几何大小相差太大，在进行流场计算时将使得 一1 2 * 西北工业大学硕士学位论文 计算过程收敛减慢，计算流场质量下降，甚至可能造成发散。由于在网格 生成过程中并没有考虑邻居网格的大小，所以在开始进行流场计算前需要 对网格的光滑度进行检查。本文采用的判断方法是要求相邻网格的加密层 之差不能大于l ；如果有，则需要将较大的网格继续拆分。 ( 2 ) 边界网格的有效性检验 一个边界网格可能与多个物面相交，如果网格内的不同物面面元之间法 向量的夹角相差太大( 超过设定值) ，此网格内在不同物面附近的流场可能 具有截然不同的流动特性，如果不把网格拆分，至少会造成流场计算的准 确度下降，还有可能导致流场计算的收敛速度减慢，甚至可能导致发散。 如图2 4 所示情形，需要把网格重新拆分。这种物面法向量判断只是一个初 步的判断，在计算了边界网格的具体边界信息后，需要对边界网格的有效 性( 关于体积) 进行判断，具体判断方法见第三章边界信息的计算。 图2 4 不规则的网格情形( 需重新拆分) 西北工业大学硕士学位论文 第三章物面边界处理与校验 直角切割网格生成之后，流场中的网格与物面将有三种不同的位置关系， 即有完整的模型内部网格、与物面相交的网格和完全在流场中的网格。如图3 1 所示。本文把完全在物体内部的网格称为固体网格，与物面相交的网格称为边 界网格，其余的称为流场网格。2 1 图3 1 直角切割网格技术中的网格类型 在直角切割网格中，边界网格的网格面与物面吻合的并不是很好。边界网 格被物面切割为流场内和物面内两部分。所谓的边界网格是指网格流场内的部 分。这样，在边界网格内的物面就变成了网格界面的一部分，而这部分网格界 面通常都是不规则的，这就造成了物面边界条件处理的困难。 考虑到离散的物面片只是计算模型的近似表示，真实的物面在物面片的交 界处也是光滑的，我们只需要求出边界网格内所有物面片的平均法向量来代表 真实物面的法线方向。有了物面的平均法向量，我们就可以完成物面边界条件 的处理。本章的主要内容是根据边界网格与物面相交的具体情况，求出边界网 格各面的面积( 有流体通过的面积) 、边界网格的体积( 其中的流体体积) 、边 西北工业大学硕士学位论文 界网格的重心和网格内物面的面积( 网格界面与物面重合部分的面积) 。另外， 本章还建立了寻找固体网格的方法，同时给出了对网格更多的有效性检验以及 对边界处理中的计算结果进行校验。 3 1 边界信息的计算 为了求得准确的流场解，就要较好的满足物面边界条件，这样要求有准确 的边界网格的几何信息。为此，当用网格切割物面后，首先求与本网格相交的 每一物面面元与网格的每一个面的交线，交线把物面分为网格内外两部分，把 在网格内的那部分保留。这样，被切后留下的物面面元都在网格内部。分别对 这些物面进行面积分，可以得到网格内每个物面的面积，它们之和就是网格内 物面的面积。以每个面元的面积作为权值，对网格内每个面元的法向量加权平 均，就可以得到网格内物面的平均法向量；在每个物面面元外法线方向且在物 面和原始网格界面之间的那部分是流场部分，在这一部分流场中分别对1 和( x ， y ，z ) 进行体积分就可以求出网格体积和流体重心。 把网格内每个物面进行拼接，取物面边缘，最后就可以得到网格面与多边 形物面的交线。这些交线可以连接成一个或者多个空间多边形，所有的边都在 物体表面。如果空间多边形有多个，则边界条件很难满足，此网格则需要重新 拆分。有一点需要强调，就是在计算交线时，线段的储存一定要有顺序，也可 以称作射线。射线的方向到物面的外法线方向在线段所在物面的投影的角为9 0 。，这样就保证原物面的外法线方向的区域都在射线的左侧，也就是流场在射 线的左侧。要计算每个面的面积，只需要把每个面上的射线按顺序连成环或者 连到边界上，这样面积就很清楚了。如图2 4 下右图所示，两匿中间部分是物 体，上下都是流场，每个面的面积通过积分就能求得。 嚣托羔啦大学硕士学证论文 3 2 删除固体潮格 在进行流场计算时，黉求将网格内的固体网格删除，这样要求对固体网格 送蟹判錾。本文爨这群判叛嚣体越揍鲍。如果边器网捂懿禁令露突全在物体爽 部，就把这个面的面积置零。这样导致的结果是，本来相邻网格间相邻西的面 积稳凌瑗在不霹了。于是，骰设粼捂a 是逑秀瓣擦，a 耧稷邻耀撩b 熬撵邻瑟 完全在物体内部，则a 的这个面丽积置为零。如果b 也为边界网格，则b 的 这个面西稷同样为零：僵翔采b 举为透秀阏捂，黉翦b 就为溺镩阚穆。不难理解， b 的非边界网格邻居c 也是固体网格，c 的非边界网格邻居d 也是固体阏格； 按此方法能够找到b 所在豹这个单连通域内的所肖固体网格。 3 3 边赛鼷格煞有效谶检验( 体积稔验) 在流场计葬串，对予个鄹俸获，j 、豹鬻格，不德会造娥特剐，l 、豹时阕多长， 而且可能在计算过程中造成奇异，减慢收敛甚至使计算发散。需毅把这种网格 并入褶邻的大网格。 本程廖采用鸵方法是把髂积较夺的嬲揍寒一个相邻懿丈瞬格关联，农计算 时把两个网格的残值求和，然后按体积平均给两个网格。具体公式如下： 耻熹r i + r j ) 巧= 熹峨鹕) 3 4 边界潮捂的封闭经校验 理论上每一个网格都楚一个封闭的嚣域，毽由于计算机精度的限制等原因 西北工业大学硕士学位论文 有些边界网格的面积在计算过程中可能出现误差，造成网格的不封闭。不封闭 的网格在进行流场计算时会造成通量项不平衡，将导致流场计算发散很快。因 此，需要对网格的封闭性进行校验。 将每个晃面的面积作为权值，对网格的所有界面的求平均法向量。理论上 封闭网格所有界面的平均法向量是一个三维的零向量。本文在计算边界网格信 息之后就给出了这个向量，这个向量就是计算误差，如果它的模较大，说明生 成的网格质量不高，可以更改生成网格的参数重新生成网格。 3 5 其它校验 为了获得高品质的网格，本文对网格还进行了以下三方面的校验。 ( 1 ) 计算模型的封闭性校验。此项校验的目的是查看所输入模型是否封闭。 与网格的封闭性校验方法相同，把每个物面片的面积作为权值，求所有面片的 平均法向量。如果平均法向量的模小于设定值，则表示模型封闭。 ( 2 ) 体积校验。此项校验的目的是查看边界网格的体积计算以及固体网格 的判断是否准确。流体的体积加上计算模型的体积是否等于流场总体积。 ( 3 ) 网格相邻面面积校验。此项校验的目的是查看边界网格的面积计算是 否准确。相邻网格相邻面的面积在计算时是分别计算的，它们是否相等也是判 断计算误差的一个方面。 西北工业大学硕士学位论文 4 1 欧拉方程 第四章欧拉方程与空间离散 本文以三维非定常欧拉方程作为求解流场的基本方程。三维欧拉方程的积 分形式5 1 如下： p d f 2 + p ( 旷a ) d s = 0 肚姬+ 肛( 旷a ) d s = - f p ( 7 , a ) d s , o v d r 2 + ( 旷a ) d s = 一j p c , a ) d s p w d f 2 + l 掣妒确d s = 一l 礤z a ) d s p e d r 2 + ( 矿a ) d s = 0 ( 4 1 ) 其中，p 、肛“、v 、w 、e 和何分别是压强、密度、直角坐标系的三个速度 分量、总能和总焓，矿= “t + v 弓+ w 乏是合速度矢量，乏，0 ，t 是沿直角坐标 轴的单位矢量。 三维欧拉方程包含五个独立的方程，但有七个未知量，方程系统不封闭。 实际求解时，必须考虑气体的热力学性质，补充两个方程，一般形式是 e = e ( p ，p )( 4 2 ) = h ( p ，p ) 对于完全气体的流动，满足下列关系式： p 砘_ 1 ) 胛一生# ) ( 4 3 ) ( 4 5 ) rrq，n r。jn r。-n rjn a一西a一西a一西a一西a一出 西乾i 照大学硕士学位论文 符= e + p p y 怒比热比，一般壤况下取1 4 。 ( 4 ，1 ) 式还可以统一地髯成下列形式： l 靛f i o d c 2 七孽刑s = o 其中， u * p 1 p u 1 一 l ，f = p , w 肛j p v 神七蟊l 珊y + p i , 掣v 专蟊z p h v 4 2 无量纲仡参数的弓i 入 ( 碡鸯 ( 4 4 ) 数值计算中，物理量一般都要无量纲化。令凡、凡、正k 分别袭示自由来流 静压强、镢度帮憨络，令互表示澈殇孛秘体鹣特徽长度。各无量缁参数定义翔 石。户，几，；= p i p 。，( i - ，i ) = ，v 。w ) ，；鬲 - g = p 。唾一h 0 | p 一爱= p 。t 毽一h 0 p 。 ( i - 三) = ( x ，y ，z ) l l ，i = f 【拓了i 】 致控方程和( 4 6 ) 式串豹物理薰替换为涎羹缡稼静携理筵，形式傈持不变。 暇翳楚厩铲瓣+ 警一学)p 。 z 英串，警= 寿+ j l 烨。2 ，坂是叁崮来流萼赫数e 西北工监大学碗士学位论文 4 3 有限体积法 本文聚用j a m e s o n 等人提出的有限体积法6 1 ，崧生成直角切割网格的 蒸礁上建立求鳃欧控方毯豹计算方法。必戴，蓥先憋( 4 4 ) 式轰接应爰于每一个 网格单元，可变形如下 翌+ 三声西西：0 西 q 未 对空润一个疆都平行予坐标辘的长方棒弼格来说，上式可离散为 詈+ 去陋+ 一毒拇；+ 蛾+ 一彭b + 幢。一霞：过b 蛰 ( 4 固 荚中，s ：、s ，、s z 分别为为长方体x 、y 、z 方囱的侧嚣积， h = p w p w u o j a w 十p p w h 下标x + 、x 一、y + 、y 一、z + 、z - 分别表示在长方体六个侧砸处的物理量。 4 毒空间离教、时闻攘进与人工耗教 ( 1 ) 空惩甏教 本文空间离散采用的中心离敞格式。( 4 7 ) 式中丘+ ，足。，o ，+ ，0 ，盈，童一都 采用中心格式进行插值。 假设瓣接i 和j 攘邻，晨网穰i 在麴格j 懿x 受辘方囊，则 髟= 露= 忙i + 声2 ，如果相邻网格大小不同，则需要按与边界面距离的大小 加权平均。如网格i 中心与边界距离为a ，网格j 中心与边界距离为b ，则 霉= 枣= 妒6 + 乒嚣k + 妨。 、il引_f 帅妣峨螂 = 一g ，；i?；j 矽 2 v w 胃 m 趣铆铆融 = 一f 西i l i 业大学硕士学位论文 这种中心差分法具有二阶空间精度。 ( 2 ) 人工耗散 使用中- t h , 格式的有限体积法，在用差分项代替微分项时，舍弃了高阶导数 项，导致某些高频误差分量在求解流场的过程中不能衰减。如果在通量项中适 当地加入一些高阶导数项，就会避免这种振荡现象，使数值解完全地收敛到定 常状态。这些人为加入的高阶项，称为人工耗散。“2 川7 1 基本耗散项是二阶差分，其系数与压强二阶差分的绝对值成正比。在激波 或驻点附近，压强梯度值很大，因此该系数也很大，而在其它区域系数则很小。 虽然二阶耗散项足以消除激波附近解的波动。但并不能使计算完全收敛到定常 状态。如果在流场的光滑部分增加四阶差分耗散项，计算可以完全收敛到定常 状态，但在激波附近加入四阶差分耗散项贝l j 会产生突跃。为了制止这种现象， 当二阶耗散项起作用时就把四阶耗散项关闭。 加入人工耗散项之后的控制方程可写为如下形式： 詈叫? 。， 肿) 一扣( o ) - d ( o ) ) 其中q ( o ) ：厅亓出，为通量项，d 为人工耗散项，r 表示残值。 赢 以历表示网格i 的的人工耗散项，d ，( 2 表示网格i 的二阶人工耗散，d f ( 4 ) 表 示网格i 的四阶人工耗散，则d 产d 产) + d y 。 研2 ) - 帮 掣) _ 秽 西2 l 。2 ) o r 。缸一谚) 弼4 ) = 一出q ，( 彰一谚) 式中，彤、彤分别为二阶、四阶自适应因子，q ，为尺度因子。 西北工业大学硕士学位论文 e = 慨一谚) ， 气，= i 蛾+ v s ，+ 船：i + c 雁丽 其中，“，v ，w ，c 分别为速度分量和当地音速，取相邻网格单元的平均值。定 义感受因子”=，则 四= k 啦m a x ( v j ，o ) s 、1 4 ，7 = m a ) 【( 0 ，k 4 1 一矗? ) 其中，弘孙，为经验常数，一般1 4 k ( 2 ) 1 ，1 2 5 6 k ( 4 ) 1 3 2 。 总的人工耗散项可表达为 d = q ，i 毋乜一谚) 一彤( 玩一玩) 】 本来是二阶中心差分，但加入人工耗散项之后，在激波附近二阶人工耗散 项起作用，精度就便成了一阶；在其他区域精度不变。 ( 3 ) 时间推进方案”1 时间推进有显式格式和隐式格式。由于隐式格式是绝对稳定格式，因此在时 间推进中没有时间步长的限制，但求解过程t e - g 复杂，耗用时间也很长。显式 格式则因为格式稳定性的限制使得时间推进中时间步长受到限制，但求解简单。 本文采用显式时间推进。 对( 4 8 ) 式，采用经典四阶龙格库塔法求解，可得以下格式 0 ”1 = 0 ”+ o l + 2 k 2 + 2 k 3 + k 4 ) 6 女t = r ( o “) k 2 - 告r ( o “+ 2 ) 岛= 告肿n + k 2 2 ) 屯= 告肿n + 屯) 西北工业大学硕士学位论文 此聒甄具硐必惭r 町1 日j 稽发，但计算近栏复杂，且比孜占用内存。米用p 咖的 简化的四阶龙格库塔格式，具有二阶时间精度，足以满足要求，计算过程简 单并能够节省内存。格式如下 0 c o ) ：0 。 0 ( 1 ) = 0 ( 。h 4 f - - i - 2 - r ( 0 ( 。) o ( 2 ) ：0 ( o ) + 竺r 巾( 1 ) 3 q 、 7 o ( 3 ) ：0 ( o ) + 竺r ( o 2 ) 2 q 、 。 0 ( 4 ) _ 俨) + 告邶) 0 ：o h 诵常只存第一步计簋人工耗散项。后面都取第一步的信。 西北工业大学硕士学位论文 第五章边界条件 求解欧拉方程，边界条件的处理至关重要。处理不当会降低计算精度，或 者计算不收敛。本文要处理的边界条件有物面边界条件、远场边界条件、对称 面边晃条件。晌1 2 1 5 1 物面边界条件 由于直角切割网格不是贴体网格，物面边界的处理极大地影响计算的精度。 欧拉方程假定流动无粘，物面边界条件就是无穿透条件，也就是物面上流动的 法向速度为零。 这样，物面上的通量 f s = 0 t p z r p z y p z ： 0 s = 0 p s , p s ， p sz 0 p 是物面上的压强，s 。，s y ，s ：是物面在三个方向的分量。 5 2 远场边界条件 计算流场时，我们要建立的计算区域不可能无限大，要准确模拟真实流动， 就要保证在计算区域远场边界处扰动不会被反射回流场。如果扰动从外边界反 射回流场，会影响收敛速度，计算精度也会大大降低。因而需要加入人为的远 场边界条件来吸收远场边界上的溢出波。本文采用一维流动的黎曼不变量来建 立远场边界条件。 西北工业大学硕士学位论文 在流动边界的法线方向存在维流动的黎曼不变量，定义如下： 霞：吼一兰 ，1 晨+ ：吼+ 墨 ，一1 瓢凑、蹑网遐发，a 刀黄愿，r 刀既然瞄，一救疆抛p y 。1 4 e 对于翘音i k , x 流，r + ，r 一均脊无穷来流计算， r + = 巩一万2 a = q 。一万2 a * 肚吼+ 鲁= 一鲁 对予怒音速爨滚，r + ，r 一均鸯诗冀城乡 援褥鬟， 晨”= 蟊一万2 a = 。一万2 a , 尺+ = 吼+ 万2 a = 一万2 a , 对于弧音速流动，r + 有计算域外插得到，r 有无穷来流计算， r = 吼万2 a = v 。一万2 a 跫+ = 蟊+ 万2 a = 一丙2 a 在入滚边暴主，边雾轰懿甥囱速度魏熵篷取熬枣来滚僮；在漤滋迭癸上， 边界点的切向速度和熵假：f i - 计算域外描得到。肖式( 4 5 ) 、( 4 6 ) 和等熵关系式 么，= 麟就磷潋诗算趱掰有鹘物瑾量。 5 3 对称瑟边赛条佟 只麓玫角没祷翻淆角的流动，所有的物理登都是关予对称瑟辩称，我稍可 以只对半模进行计算，不但可以带省内存，而且节省了一半的计算时间。对称 2 5 * 西北工业大学硕士学位论文 面边界可以使用镜面反射边界条件，横跨对称面的变量都可以简单的通过对称 面镜面反射求得。本文采用的就是镜面反射对称面边界条件。 西北：泣大学硬士学位论文 第六章加速收敛技术 采潮雅定常敬拉方箍和裙篷条 牟求瓣定常弱鞭时，不要求完企模叛真实懿 非定常过程，可