（流体力学专业论文）非结构多重网格法流场数值模拟.pdf

南京航空航天大学硕士学位论文 摘要 l 在流体力学数值模拟的领域中，多重网格法被公认为一种有效的加速收敛措施。 l 本文对非结构多重网格法进行了研究，其中包括重叠式非结构多重网格法和聚合式非 结构多重网格法，并且分别给出了相应的算例。特别对于聚合式非结构多重网格法， 本文设计了一种新的非结构网格的聚合方法，使得聚合成的粗网格具有很好的质量。 在多重网格法求解e l l l c r 方程的时候，在每一重网格上迭代的时候我们采用有限体积 方法。 ， f 此外，针对本文多重网格法方面的研究，在本文最后介绍了相关的软件开发及数 据可视化方面的工作。4 关键词：多重网格法；非结构网格；有限体积法i 欧拉方程；软件设计x 数据可视化a 非结构多重网格法流场数值模拟 a b s t r a c t i nt h ef i e l do fn u m e r i c a ls i m u l a t i o no ff l u i df l o w , m u l t i g r i da l g o r i t h mi sae f f e c t i v e w a y t oa c c e l e r a t et h ec o n v e r g e n c e i nt h i s p a p e r , a p p l i c a t i o n so ft h ep r i n c i p l e so ft h e m u l t i 【g r i da p p r o a c ht o u n s t r u c t u r e dm e s hp r o b l e m sa r ed e s c r i b e d b o t ht h eo v e r l a p p i n g m u l t i 【g d da l g o r i t h ma n dt h ea g g l o m e r a t i o nm u l t i g r i da l g o r i t h mh a sb e e ni l l u s t r a t e da n d e x a m p l e so fp r a c t i c a la p p l i c a t i o n s h a v eb e e ng i v e n e s p e c i a l l yf o rt h e a g g l o m e r a t i o n m u l t i g r i da l g o r i t h m ，an e w r o b u s tm e t h o dt og e n e r a t ec o a r s eg r i d sh a sb e e n d e v e l o p e d ，b y w h i c hw ec a ng e n e r a t ec o a r s em e s h e sw i t hg o o dq u a l i t y w h e na p p l y i n gt h em u l t i g r i d a l g o r i t h m ，t h ef i n i t ev o l u m e m e t h o di sa d o p t e do ne a c h l a y e r t os o l u t et h ee u l e r e q u a t i o n s b e s i d e s ，i ns u b s e q u e n ts e c t i o n s ，s o m e r e s e a r c h e sa b o u tt h er e l a t e ds o r w a r e e n g i n e e r i n g a n dv i s u a l i z a t i o no f d a t aa r ed e a l tw i t h k e y w o r d s ：m u l t i g r i da l g o r i t h m ；u n s t r u c t u r e dg r i d ；f i n i t e v o l u m e m e t h o d ；e u l e r e q u a t i o n s ；s o f t w a r ee n g i n e e r i n g ；v i s u a l i z a t i o no f d a t a ； l l 南京航空航天大学硕士学位论文 第一章前言 自上世纪后期，随着计算机技术的飞速发展，计算机已经成为人类解决问题 的重要工具。借此东风，c f d ( 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 方法以其花费少、周期短的特点成为除实验方法外的又一重要途径。 虽然c f d 技术已经成为研究流体运动的重要手段，但到目前为止它还不能 够在实际运用中完全取代实验方法，其在许多方面还有待完善。需要提高的方面 主要包括计算的精度和计算速度。在提高精度方面，国内外从事c f d 研究者先 后提出了各种高精度格式，实践证明，取得了良好的效果。在提高速度方面，国 际上先后出现了当地时间步长、残值光顺、隐式格式、多重网格法等一系列加速 收敛的措施，其中新兴的多重网格法在近几年引起了广泛的关注，虽然目前对它 的研究尚没有达到非常成熟，但是多重网格法已经被公认为非常有效的加速收敛 手段。迄今为止，国际上c f d 领域的多重网格法研究已经有了相当的进步，以 d j m a v r i p l i s 为代表的部分科技工作者已经做出了相当有说服力的算例。 c f d 领域的多重网格法主要是为了加速收敛速度，节约计算时间，并且使 得整个收敛过程比单重算法明显更加稳定。目前国际上对c f d 领域的多重网格 法的研究主要分两个分支：基于结构网格的多重网格法和基于非结构网格的多重 网格法。其中对于基于结构网格的多重网格法的研究起步较早，也相对比较成熟。 主要因为结构网格在网格的结构上层次清晰，点与点以及元素与元素之间的结构 关系非常明确且很有规律。因此在结构网格上来实现多重网格算法时，各层网格 之间的拓扑关系非常明确，粗、细网格之间的元素数目比非常容易控制，相邻两 层网格之间的数据传递关系容易实现。由于以上原因，大多数的研究者采用结构 网格上的多重网格法。但是随着c f d 技术的迅速发展以及其应用领域的不断拓 宽，非结构网格由于其结构上的灵活性而被广泛采用，特别在对复杂外形的处理 上非结构网格显示出了比结构网格更大的优越性。基于以上原因，如何在非结构 网格上实现多重网格算法成为了一个新的课题。 目前国际上对非结构多重网格法的研究分为两个分支：重叠式多重网格法和 聚合式多重网格法。但无论是重叠式还是聚合式，都面临着比结构多重网格法更 大的难度，这主要是由于非结构网格结构上的随意性所带来的粗网格生成以及 粗、细网格之间数据传递关系的确定这两方面的难度。目前在国际上，非结构多 重网格法在c f d 领域已经由了比较成功的应用，但是还不是特别成熟。在国内已 经有人做了一些为数不多的尝试，但是大多是二维情况下的二重多重网格的方法 性尝试，而且效果还远没有达到要求，特别是在四重非结构多重网格方面，至今 我们还没有发现在国内有非常成功的算例。 本文对于重叠式多重网格法和聚合式多重网格法都做了相应的研究。特别对 于聚合式多重网格法，本文分别针对二维和三维的情况自行设计了新的网格聚合 方法。这种新型的聚合思路不仅聚合的粗网格质量好，而且聚合的速度非常快， 无论对于二维还是三维的情况，其聚合网格所耗费的时间相对于迭代计算的时间 非结构多重网格法流场数值模拟 来讲是微不足道的，这就保证了网格的聚合过程基本上不占据多重网格算法的时 间，保证了聚合式多重网格法相对于单重算法在收敛速度上的优越性。 另外，本文为了适应当今软件开发的潮流，对于本文所提到的多重网格法进 行了相应的软件封装，使得程序运行时的人机交互更加直观和简便。针对本文流 场数值模拟的工作，本文分别对于二维和三维情况下的流场数据进行了相应的可 视化方面的研究，使得我们可以更加直观的对各种数据进行分析。 南京航空航天大学硕士学位论文 第一章e u l e r 方程的有限体积解法 本章将介绍在非结构网格上求解e u l e r 方程的有限体积解法“”“1 。该方法 采用格心格式对方程进行空间离散，并且离散成的网格单元可以是任意形状。 2 1 控制方程 守恒型二维e u l e r 方程可以写为如下形式： 昙p + ( 肋一强 - ， 其中，x ，y 为坐标，w 为流场变量。q 为边界为s 的控制体。 f ，g 为流动通量。 f = p u p u 2 + p p u v p u h w = p p u p v 施 g = p v p u v p v 2 + p p v h ( 1 2 ) ( 1 3 ) 沿着边界的线积分是按逆时针方向进行的。 p ，p ，h 和e 分别为密度、压强、焓和能量，u 和v 分别为沿着坐标轴方向 的速度分量。对于完全气体： p e = e ( ，一1 ) + p ( u 2 4 - v 2 ) 2 ( 1 4 ) n 日= p e 4 - p 2 2 有限体积空间离散 ( 1 5 ) 有限体积的空间离散首先将计算域无重叠的单元，然后在每个单元上积分守 恒方程。因为单元不随时间变化，因此方程1 1 可以写为如下形式： 里：一；( f d y ；- g d x )( 1 ，d ) 一r 一 1 ， a tl m 。 非结构多重网格法流场数值模拟 其中q 和s 分别代表该单元的面积和边界，w 为单元上的变量的平均值。 对公式1 6 进行空间离散后可以写成如下形式： 警一q , 。n 。 ( 1 7 ) q 。为第k 个单元的面积，w - 为该单元上的守恒变量。q - 为通量积分的离散逼近， 形式如下： j g = z - ( f a y - g a , z ) ( 1 8 ) l i lt 其中 薯= 一x 。 缈。= y 6 一y 。 公式( 1 9 ) 中下标a ，b 分别表示边的两个顶点。 图1 1 相邻元素示意图 f i g 1 1 n o m e n c l a t u r ef o rc e l lc e n t e r sa n dg r i dv e r t i c e s ( 1 _ 9 ) 其中求和是在构成该单元的边上进行的。由于本文的有限体积法是格心格式的， 因此各个守恒变量是代表单元形心的值，而边上的值我们由该边左右两边的单元 上的值的平均求得，如下式： 彬= ( + ) 2 另外，在我- f l v j i 入 z f = u ，a y ，一巧缸 后，单元k 上离散形式的e u l e r 方程如下： ( 1 i 0 ) ( 1 1 1 ) 4 一塑室塾窒堕蒌奎兰堡主堂堡堡塞 盟：一土事 d t q 。智 z 白 z p u 七p 姆 z d 矿一p a x z p h 上面的一般形式即可以用在结构网格上也可以用在非结构网格上。 2 3 人工耗散 ( 1 1 2 ) 上述格心格式是无耗散的，离散误差等没有得到阻尼，因此可能出现振荡。 为了消除振荡，我们在方程的右边加入人工耗散项，对于单元k ，方程变为： 警一( 幺一o d c a 。 ( 1 1 3 ) 其中耗散项d t 由二阶和四阶两部分组成，四阶部分在解比较光滑的区域起作用， 在激波附近将被关闭；= 阶部分在激波附近将被激活来抑制振荡。其中打开和关 闭的开关是通过一压强变化的感知器来实现的。d 。的形式如下： n a g 虹d 龋5 d i = d j 2 + d j 4 ( 1 1 4 ) 其中 d j 2 = 口0 2 ( 阡0 一) d ；4 = 一口0 4 ( v 2 降0 - v 2 ) 。 ( 1 1 5 ) 上式中，i 为分割元素k 和p 的边，v2 定义如下： v 2 = ( 一职) ( 1 1 6 ) 自适应系数e 。和e “。定义如下： 占j 2 ) = | ( 2 ) y 其中 0 4 lm a x ( o ，七f 4 ) 一0 2 ) ( 1 1 7 ) v t = ( 1 1 8 ) 比例因子a 。定义如下： q = u a y 一险l + c 瓜葡 ( 1 1 9 ) 其中u ，v ，c 为边上的平均值。c 为当地音速。 最一最 一一+ 0 0 非结构多重网格法流场数值模拟 2 4 时间推进 我们通过对微分方程时间上进行积分来达到稳定状态，形式如下 娶a t 咄 上式右边代表单元k 上的残值误差。 r 。= 一( g d k ) n i 本文中我们采用显式四步格式进行时间推进： 矽蚋= w “ 其中 2 。5 边界条件 2 5 1 物面边界条件 w 伽) = w ( o ) + 口m a t r t ”一1 ) m = l t 0 4 ”1 = 矽4 1 5 i 1 口2 = ： 1 口3 = j o t 4 = 1 ( 1 2 0 ) ( 1 2 1 ) ( 1 l2 2 ) ( 1 2 3 ) 本文中我们在物面采用无穿透边界条件“”。在具体的实现上只需令公式 1 1 1 中的z i 为0 即可。 2 5 2 远场边界条件 本文在远场处采用无反射边界条件 1 9 1 。对于远场边界不能简单的规定自由 流边界条件，否则流场内部的扰动波传到远场边界会反射回来，对解的收敛不利。 沿流场边界法向我们作一维特性分析，得到r i e m a n n 不变量。现用下标o 。表 示自由流值，下标e 表示外插得到的值。又令香，a ，c 分别为速度向量、单位 法向量和音速。则得到r i e m a n n 不变量如下： 6 南京航空航天大学硕士学位论文 耻扎瓦一鲁 耻酗覆+ 鲁 由公式( 1 2 4 ) 和( 1 2 5 ) 得到： 牙蠢= 丢( r 。+ r 。) c ：r - l ( r 。一只。) 4 9 ( 1 2 4 ) ( 1 2 5 ) ( 1 2 6 ) ( 1 2 7 ) 若远场边为出流边界，则切向速度分量和熵值取外插值；若远场边界为入流 边界，则切向速度分量和熵值取来流值。对于超音速流场，若为如流边界，则 矿= k ；若为出流边界，则= 睨。 非结构多重网格法流场数值模拟 第二章重叠式非结构多重网格流场数值模拟研究 非结构网格上的多重网格法思想“”“”1 分为重叠式和聚合式两 种，重叠式非结构多重网格法“1 式在一系列粗网格和细网格之间进行值的传递， 并且粗网格和细网格结构无关，各层网格单独生成。该方法层与层之间数值传递 的难度相对较大。本章采用重叠式多重网格法对n a c a 0 0 1 2 翼型的亚、跨两个状 态都进行了计算，在每一层网格上采用有限体积法“2 “1 进行计算。通过与 单重算法的效率比较，多重网格的高效性得到了充分的体现。 2 1 采用的网格 8 南京航空航天大学硕士学位论文 图2 1 采用的非结构网格 f i g 2 i t h ea d o p t e du n s t r u c t u r e dm e s h e s 本文的重叠法0 1 采用了四重结构无关的非结构网格，见图2 1 。所谓结构无 关是指各层网格独立生成。本章中，我们对n a c a 0 0 1 2 翼型采用d e l a u n e y 方法”1 分别生成了由密到粗的四重非结构网格，这四重非结构网格除了内、外边界相同 之外，其他结构方面完全无关。 独立生成相邻的粗网格和细网格的单元数比约为1 ：4 ，最密的第一重网格的 单元数为9 3 4 0 ，结点数为4 7 6 7 。 9 非结构多重网格法流场数值模拟 2 2 数值传递 图2 2 粗网格和密网格示意图 f i g 2 2 t h eo v e r l a p p i n gm e s h e s 相邻网格层之间的参数的传递采用面积加权2 1 的思想。从细网格到粗网格 的数值传递如下： 睨= ( x s i c ) s 。 ( 2 1 ) r 。= ( r ，x s i o ) s ， ( 2 2 ) 其中w 表示流动变量，r 表示残值，s 表示单元面积，而s 。则表示粗网格单元与 细网格单元相交部分的面积。如上图所示。下标。和，代表粗网格和细网格。 粗网格到细网格的修正量插值为： c ，= ( c 。x s i c ) 丹 ( 2 3 ) c 。与c ，分别表示粗网格和细网格的修正量。 2 3 多重网格算法 本文采用的方法的基本思想是在细网格和粗网格上交替迭代“3 “”“”“”“” “”，循环方式见图2 3 。图中第一层网格为初始最密的网格，第二到第四层网格 依次逐渐变疏。 用粗网格上的迭代过程来消除细网格上的低频误差，以此来提高收敛速度。 本文中在每一重网格上都采用四步显式r u n g e - k u t t a 法“玎“1 进行时间推进。 循环具体过程如下： ( 1 ) 首先在最密的第一层网格上进行一次迭代。 1 0 南京航空航天大学硕士学位论文 w 9 = w ”一a 。值( 4 1 ) ( 2 4 ) ( 2 ) 将流动参数传递到第二层粗网格上，粗网格上进行一次迭代。此时需加一 强迫函数： p = 胄一r ( w 、 w 们= ”一f k ( 州) + p 】 ( 2 5 ) ( 2 6 ) r 和w 为从上一层网格传递过来的参数。然后依次将值传递到更高层的粗网格 上进行迭代。 ( 3 ) 最高一层粗网格迭代完后，将修正量直接回插到最密的第一层网格上，此 过程中各层网格上不再进行迭代。到此作为一个循环。 ( 4 ) 然后再返回第一步、如此反复，直至收敛。 2 4 算例分析 图2 3 循环过程 f i g 2 3m u l t i g r i dc y c l e 本文对亚音速( 来流马赫数m = 0 2 ，迎角a = 1 2 5 。) ，跨音速( 来流马赫数 m 。= o 8 ，迎角a = 1 2 5 。) 两个状态分别用多重网格及单重网格( 单独用最密一 层) 进行了计算。 图2 4 为亚音速情况下( 来流马赫数m = 0 2 ，迎角a = 1 2 5 。) 收敛史的比 较。图2 5 为跨音速状态下( 来流马赫数m = o 8 。迎角a = 1 2 5 。) 收敛史的比 较。图中“s i n g l e ”表示单重算法的收敛历史；“m u l t i ”表示多重网格算法的收 敛历史。图2 4 和图2 5 中横坐标i t e r a 表示循环数。对于多重算法，每一个循 环表示在各层网格上都进行一次r u n g e k u t t a 迭代“1 。对于单重算法，一个 循环就代表一次r u n g e k u t t a 迭代。而纵坐标代表残值误差r e s 的对数( 以1 0 为底) ，即误差的量级。残值误差的定义如下： 斐堕塑兰重旦整鲨亟堑鍪堡堡塑 r e s ； 10 0 o0 0 10 0 一二j u j0 0 40 0 i l 】c 00 0 1 一i ：0 c ? - 一。 l i 一一二o o 40 0 5o o 60 0 10 0 o 套0 0 d j 叫m r 图2 4 m = 0 2 ，= 1 2 5 。时的收敛史 f i g 2 4 t h e h i s t o r yo f c o n v e r g e n c e w h e n m = o 2 a = 1 2 5 。 图2 5m = 0 8 ，q = 1 2 5 。时的收敛史 f i g 2 5 t h e h i s t o r yo f c o n v e r g e n c e w h e n m = 0 8 a = 1 2 5 。 十er - 0 o o 南京航空航天大学硕士学位论文 由以上两个图可以看出多重网格法的高效性：m 。= o 2 时，残值达到l o - 5 量级 需3 0 0 多个循环左右；m = 0 8 时，残值达到1 0 量级需6 0 0 个循环左右。 用多重网格计算时，粗网格上的计算耗时约是其相邻密网格上耗时的i 4 ， 以此推算，各个状态下多重和单重残值降到1 0 5 量级的耗对比如下：( 1 ) m 。= o 2 时，约为l ：1 2 ( 2 ) m 。= 0 8 时，保守估计为1 ：3 5 。 由以上两个收敛史图象可以看出多重网格法的显著效果。跨音速情况下，其 抑制振荡的作用比较突出；亚音速情况下，多重网格法可以非常有效的解决单重 情况下残值停滞不前或衰减非常缓慢的问题，使得残值能很快的衰减到l o 。量级 以下。 图2 6 和图2 7 各状态下单重算法和多重算法的结果比较，由图可见，结果 完全重合。图2 6 和图2 7 中上、下两条曲线分别表示下、上翼面压强系数分布。 图2 6m = 0 2 ，a = 1 2 5 。时的c ，曲线 f i g 2 6 d i s t r i b u t i o no fc 。 w h e nm = 0 2 n = 1 2 5 。 图2 7m = 0 8 ，a = 1 2 5 。时的c p 曲线 f i g 2 7 d i s t r i b u t i o no fc - w h e nm = o 8 a = 1 2 5 。 l 蛐 。 i j o 山 一 叫 。一 非结构多重网格法流场数值模拟 2 5 结论 本章用非结构重叠式多重网格法对n a c a 0 0 1 2 翼型的亚、跨状态进行了无粘 绕流的数值计算，并与单重网格计算的效率进行了分析和比较，验证了多重网格 高效的加速收敛效果。由于用多重网格计算时，粗网格上的计算用时相对很少， 因此可以节约大量的机时。 多重网格法之所以能够如此大幅度的加快收敛速度，取决于多重网格法机制 对误差的显著消除效果。我们知道用单重算法的时候，高频误差的消除比较容易， 但是低频误差的消除则比较困难。但是当采用多重网格法的时候，初始密网格上 的低频误差在粗网格上就变成了高频误差，因此多重网格法的思想就是利用一系 列粗网格来加速初始密网格上低频误差的消除。这样在粗、细网格上交替迭代就 能够使得初始密网格上的高频误差和低频误差都得到很快的消除。因此多重网格 法比单重算法体现出了非常明显的收敛速度和收敛稳定性上的优势。 1 4 南京航空航天大学硕士学位论文 第三章二维非结构聚合式多重网格法数值模拟研究 对于聚合式非结构多重网格法“1 ，其思想是由初始的非结构网格聚合成几套 形状不定的多边形粗网格，在生成粗网格的过程中自然地记录各层之间的传递关 系。聚合韵方法决定了网格的质量和各层之间的网格数目比。聚合法的优点是各 层网格之间的数值传递关系较重叠法简单，但如何聚合高质量的粗网格是问题的 关键。本章设计了一种高效的网格聚合方法，该方法的聚合过程采用层推进的思 想，聚合过程完全自动化，并且能够很好的保证相邻粗、细网格之间的网格单元 数目的比例关系，从而保证了多重网格算法对误差的消除效果。本章采用四重聚 合式多重网格法分别对n a c a 0 0 1 2 翼型和三段翼型进行了二维e u l e r 方程的数值 求解，并与相同c f l 数下的单重网格的计算效率进行了比较，体现了多重网格法 加速收敛的高效性。 3 1 网格的聚合生成 在上一章介绍的重叠式多重网格中，采用的各层由密到疏的网格是单独生成 的，由于每一层网格都是三角形的非结构网格，因此从单重网格的质量上来讲比 较容易保证。而作为聚合式多重网格法，粗网格是由密网格聚合生成的。从形状 上来讲，粗网格单元是不规则的多边形，网格质量的优劣完全取决于聚合的方法。 因此如何快速的聚合生成一系列高质量的粗网格将是本章的一项重要内容。 作为聚合式多重网格技术的关键部分，网格聚合质量的优劣无疑对计算效果 有着举足轻重的影响。因为多重网格法的思想就是在密网格上消除高频误差，而 密网格上的低频误差由于很难消除，因此采用更粗的网格来帮助消除密网格上的 低频误差。因此粗细网格之间的尺度关系对于误差的消除具有决定性的影响。本 章设计了一种高效的网格聚合方法，由该方法聚合生成的粗网格质量好，而且聚 合过程的耗时很少。用奔腾1 1 1 5 5 0 的p c 机聚合如图1 所示的四套多边形网格只 需3 秒左右( 包括读写文件的时间) 。其中初始的n a c a 0 0 1 2 翼型三角形网格的网 格数为9 3 4 0 ，节点数为4 7 6 7 ，多段翼型网格数为1 2 1 3 8 ，节点数为6 2 7 3 。各层 相邻粗网格的之间的网格数比都约为4 ：l 。本文中的网格聚合过程如下： ( 1 ) 首先以每个节点为中心，以此节点的相邻三角形单元的形心和边的中 点的连线围成的封闭多边形作为第一层粗网格的单元。因此第一层粗 网格单元数与初始网格节点数相同。 ( 2 ) 由第一层粗网格生成第二层粗网格采用了一种类似于层推进的方法： 首先，选择内外边界上的元素作为起始聚合边界，并在聚合边界上任 意选择一元素作为待聚合元素。将与此待聚合元素相邻的没有聚合过 的密网格的元素全部聚合为一个粗网格，同时将所有与刚聚合过的密 网元素相邻的未聚合的密网元素记录到待聚合边界中。然后选择聚合 边界中与刚生成的这个粗网元素相邻的密网元素作为待聚合元素，将 与此待聚合元素相邻的没有聚合过的密网格的元素全部聚合为一个 粗网格如此反复，直到聚合边界中所有的密网元素都聚合完毕。 非结构多重网格法流场数值模拟 然后将聚合边界更新为待聚合边界，选取待聚合元素重复以上步 骤直到所有的密网元素都聚合完毕。此时，第二层粗网格生成完毕。 ( 3 ) 其它层粗网格的生成与第二层粗网格的生成过程完全相同。 本章采用上述方法对n a c a 0 0 1 2 翼型进行了粗网格的聚合生成，结果如 图3 1 。 6 南京航空航天大学硕士学位论文 非结构多重网格法流场数值模拟 ( e ) 图3 1 本文中采用的n a c a 0 0 1 2 翼型网格 t h ea d o p t e du n s t r u c t u r e dm e s h e sa r o u m dn a c a 0 0 1 2a i r f o i l 为了验证本文设计的网格聚合方法的通用性和鲁棒性，本文还对形状相对更 加复杂的三段翼型进行了租网格的聚合生成，结果如图3 , 2 。 南京航空航天大学硕士学位论文 非结构多重网格法流场数值模拟 ( e ) 图3 2 本文中采用的多段翼型网格 f i g 3 2 t h e a d o p t e d u n s t r u c t u r e dm e s h e sa r o u n dt h em u l t i e l e m e n ta i r f o i l 3 2 数值传递 相邻网格层之间的参数的传递采用面积加权“”儿8 ”“ 细网格到粗网格的数值传递如下： 形= ( x s 扣) 瓜 2 1 的思想。从 ( 3 1 ) r 。= ( r ，s 扣) s ， ( 3 2 ) 其中w 表示流动变量，r 表示残值，s 表示单元面积，而s 。则表示粗网格单元与 细网格单元相交部分的面积。如上图所示。下标。和，代表粗网格和细网格。 粗网格到细网格的修正量插值为： c i = ( c o ) s ， ( 3 3 ) c 。与c ，分别表示粗网格和细网格的修正量。 需要说明的是，以上公式中的s ，。实际上就是细网格单元的面积。因为粗网 格是由细网格聚合生成的，细网格单元完全包含在它所从属的粗网格单元的内 部。 2 0 南京航空航天大学硕士学位论文 图3 3 粗、细网格的结构关系 f i g 3 3 t h ec o a r s em e s ha n dt h ef i n em e s h 由于本文采用聚合式粗网格，如图3 3 ，因此粗、细网格之间的结构关系比 较明确，数据传递的实现也相对于重叠式多重网格法更加容易。如图3 3 所示， 每个密网格单元都完全包含在它所从属的粗网格单元中，这种情况下，粗、细网 格单元相交部分的面积就是密网格单元的面积。这样粗、细网格之间的数据传递 关系就更加明确，而不必像重叠式多重网格法一样还要花费时间求粗、细网格单 元相交部分的面积。这是聚合式非结构多重网格法的优点所在。 3 3 循环迭代 在以往的多重网格法的研究中，人们大多采用v 型循环“1 和w 型循环“，。 一般来讲。v 型循环过程比较简单，w 型循环相对比较复杂。一般情况下w 型循 环达到收敛所需的循环数要少，但同时由于其每一个循环在各层粗网格上需要更 多次的迭代，因此也需要在反复的插值过程中浪费时间。因此可以说，这两种循 环在都能达到较好效果的基础上各有利弊，对于不同的算例。很难断定哪一种方 法更优越。一般来讲，采用w 型循环是为了更彻底地在粗网格上消除误差，但一 般情况下，v 型循环已经能够取得较好地效果，因此目前的多重网格法的应用主 要是采用v 型循环，当然v 循环达到收敛所需要的循环数要稍多一些。本文考虑 到以上原因，采用了一种折中的方案。循环过程如图3 4 所示。考虑到最粗的两 层网格的网格数和边数都很少，因此本文只在最粗的两层网格上重复迭代，这样 每一个循环比w 型的一个循环节约很多机时，另外，达到收敛所需的循环数要比 v 型循环少。 图3 4 循环过程 f i g 3 4m 蜊do y e l e 在最密的第一层网格上进行迭代公式如下： w 9 = 们一口。吠( 形9 1 ) 在粗网格上进行迭代时公式如下： p = r 一r ( w ) ( 们= 矽( o ) 一f k ( 形( 州) + p 】 ( 3 4 ) ( 3 5 ) ( 3 6 ) r 和w 为从上一层网格传递过来的参数。然后依次将值传递到更高层的粗网格 上进行迭代。 最高一层粗网格迭代完后，将修正量直接回插到最密的第一层网格上，此过 程中各层网格上不再进行迭代。到此作为一个循环。然后再返回第一步o 如此反复，直至收敛。 3 。4 算例分析 本文对跨音速( m 。= o 8 ，a = 1 2 5 。) 状态下n a c a 0 0 1 2 翼型和亚音速( m 。= 0 3 ， a = 1 0 。) 状态下的多段翼型的绕流分别进行了单重算法和多重算法的数值模拟。 本文采用的网格是第一到第四层粗网格( 图3 1 和图3 2 中的b ，c ，d ，e ) ， 四层粗网格的网格数基本上都是按4 ：1 递减的。 本文对单重算法及上文提到的介于v 型与w 型之间的折中多重算法的收敛速 度进行了比较，结果如图3 5 和图3 6 ，其中横坐标表示循环数，纵坐标表示残 值误差的对数，即收敛的量级。图3 7 和图3 8 是两种翼型的多重网格算法和单 重算法计算的c ，曲线对比，由图可见，结果完全重合。 为了便于比较，本文中单重算法和多重算法在在最密一层网格( 即图3 1 和 图3 2 中的b 网格) 上取相同的c f l 数。对跨音速( m = 0 8 ，= 1 2 5 。) 状态 下n a c a 0 0 1 2 翼型用多重算法只需3 2 0 个循环便可以收敛到1 0 5 量级，而单重算 法需2 0 0 0 步以上。 为了进一步验证这种多重网格算法的高效性本文还对多段翼型的绕流进行 了数值模拟。这种多段翼型由于其形状的原因，通常采用显式单重算法收敛很慢。 但亚音速( m = 0 ，3 ，a = 1 0 。) 状态下的多段翼型用折中算法只需8 0 0 个循环左 右，单重算法需5 4 0 0 步以上。 一壹矍堕窒堕丕盔兰堡主兰垡堡奎 40 0 00 0 i 】 l i 止一4 0 0 一 i 二r l 一一80 0 一 一 二0 0 一1 60 0 一二00 0 0 00 0 l 。 u _ _ _ ，一二o 。 二 j 4o o 图3 5 收敛历史( n a c a 0 0 1 2 翼型) f i g 3 5c o n v e r g e n c e r a t eo f a g g l o m e r a t i o n m u l t i 鲥d a l g o r i t h mc o m p a r e d w i t h s i n # - g d da l g o r i t h m ( n a c a 0 0 1 2a i 响i l ) 图3 6 收敛历史( 多段翼型) f i g 3 6c o n v e r g e n c e r a t eo f a g g l o m e r a t i o n m u l t i g r i d a l g o r i t h mc o m p a r e d w i t hs i n g l - g r i da l g o r r h m ( m u l t i - e l e m e n ta i r f o i l ) 1 1 _ 三r 、0 0 0 斐堕塑墨蔓塑整鲨煎塑墼堡堡垫 o o ) 5 0 o0 0 o5 0 0 0 6 0 03 0 04 0 o4 0 一o8 0 图3 7c 。分布( n a c a 翼型) f i g 3 7 d i s t r i b u t i o no fg w h e nm 。= 0 8 a = 1 2 5 。 ( n a c a 0 0 1 2a i f f o i l ) 图3 8c 。分布( 多段翼型) f i g 3 8 d i s t r i b u t i o no fc ， w h e nm - = 0 3 ，a = l o 0 。 ( m u l t i - e l e m e n ta i r f o i l ) 南京航空航天大学硕士学位论文 3 5 结论 本文设计并采用了一种新的网格聚合方法，此方法的基本愚想是以单元为基 本元素组成待聚合边界进行推进式的聚合。此方法的优点在于能够保证相邻两层 网格的单元数比约为4 ：1 ，并且聚合的网格具有较好的质量。另外，本文采用了 一种介于v 型循环和w 型循环之间的循环方式。此方法只在最粗的两层网格上重 复迭代，由于这两层边数和网格数都很少，额外的迭代并不花费多少时间，因此 每个循环的耗时比w 型少的多，但又比单纯的v 型循环较明显的节省循环数。通 过以上算例，本文验证了此种多重网格法加速收敛的高效性。 非结构多重网格法流场数值模拟 第四章三维非结构聚合式多重网格法数值模拟 对于三维聚合式非结构多重网格法口舶1 1 ”，由于非结构网格结构上的任意性， 对网格聚合的控制需要更多的技巧。本章将多重网格法推广到三维非结构网格的 情况，自行设计了一种高效、高质量的三维非结构网格聚合方法，并采用四重非 结构网格进行三维e u l e r 方程聚合式多重网格法的数值求解。通过与单重网格的 计算效率进行的比较可以看出本文将多重网格法运用到三维非结构网格的尝试 是非常成功的。 4 1 三维非结构网格的聚合 本文白行设计了一种高效的三维网格聚合方法，在保证网格质量的同时，也 充分考虑到了聚合速度和聚合算法的鲁棒性。 本文中初始四面体网格的单元数为2 9 1 9 1 4 ，在此基础上聚合三重粗网格所 需要的时间不到3 0 秒( 2 5 6 兆内存的3 5 0 p c 计算机) 。相邻的两重粗细网格的网 格数比约控制在1 ：8 左右，三重粗网格的网格数分别为5 0 3 2 5 ，5 7 1 l ，7 9 7 。 具体的聚合方法如下： ( 1 )选内外边界的点组成待聚合点集合，然后选其中一点为聚合点，将所有 包含此点的元素聚合为一个粗元素。然后从待聚合点集合中选取一个与原 起始点相邻的点为聚合点，将所有包含此点的元素聚合为一个粗元素。如 此反复，直至待聚合点集合中的所有点都被聚合完为止。用位于流场内部 的新生成的一层粗网格表面的点来更新待聚合点集合，然后选一点为聚合 点如此反复，直到聚合完毕。这样就由第一层细网格生成了第二层粗 网格。 ( 2 )将第二层粗网格中的元素两两聚合为一个元素，这样聚合成一层过渡网 格。 ( 3 )选过渡网格中边界上的元素组成起始聚合边界，并在聚合边界上任意选 择一元素作为待聚合元素。将与此待聚合元素相邻的没有聚合过的密网格 的元素全部聚合为一个粗网格，同时将所有与刚聚合过的密网元素相邻的 未聚合的密网元素记录到待聚合边界中。然后选择聚合边界中与刚生成的 这个粗网元素相邻的密网元素作为待聚合元素，将与此待聚合元素相邻的 没有聚合过的密网格的元素全部聚合为一个粗网格如此反复，直到聚 合边界中所有的密网元素都聚合完毕。然后将聚合边界更新为待聚合边界， 选取待聚合元素重复以上步骤直到所有的密网元素都聚合完毕。此时， 第三层粗网格生成完毕。 ( 4 )重复( 2 ) ，( 3 ) 步骤，生成第四层粗网格。 从上述的聚合步骤可以看出，本文提出的这种聚合思想具有较高的智能化和 自适应特性。内边界和外边界形状上的复杂度对聚合过程并无影响，也就是说， 本文所设计的三维网格聚合方法具有相当好的通用性，聚合过程中不会出现类似 于单重网格生成中的网格生成失败的情况。因此本文设计的三维网格聚合方法本 南京航空航天大学硕士学位论文 身就是针对于任意的三维复杂外形。 4 2 数值传递 相邻网格层之间的参数的传递采用面积加权“”们”“”1 的思想。从 细网格到粗网格的数值传递如下： 睨= ( x s ，c ) s 。 ( 4 1 ) r 。= ( r ，s 弦) s ， ( 4 2 ) 其中w 表示流动变量，r 表示残值，s 表示单元面积，而s 。则表示粗网格单元与 细网格单元相交部分的面积。如上图所示。下标。和，代表粗网格和细网格。 粗网格到细网格的修正量插值为： c ，= ( c o s 斥) 彤 ( 4 3 ) c 。与c ，分别表示粗网格和细网格的修正量。 与二维聚合式多重网格法相同，以上公式中的s 。实际上就是细网格单元的 面积。因为粗网格是由细网格聚合生成的，细网格单元完全包含在它所从属的粗 网格单元的内部。 4 3 循环迭代 本文采用v 型循环n 。首先在最细的第一层网格上进行一次显式r u n g e k u t t a 迭代“1 ，然后将数据传递到第二层网格，在第二层网格上进行一次 r u n g e k u t t a 迭代直至第四层。第四层上迭代完毕后，将修正量直接传 递到第一层网格。如此反复，直至收敛。 图4 1 循环过程 f i g 4 1m u l t i g r i dc y c l e 非结构多重网格法流场数值模拟 本文中在每一重网格上都采用四步显式r u n g e - k u t t a 法l “i 进行时间推进。 在最密的第一层网格上进行的迭代为： 形。= 矽”一口。氓( 形9 1 ) 租网格上进行迭代时，需加一强迫函数以修正传递误差： p = r 一r ( 矿) ( 们= ( ”一a q a t r ( w g 一1 ) + p 】 ( 4 4 ) ( 4 5 ) ( 4 6 ) r 和w 为从上一层网格传递过来的参数。式( 5 ) 中r 和w 分别表示由密网格 传递过来的残值和流动变量。 4 4 算例分析 本文对跨音速( m 。= o 8 4 ，a = 3 0 6 。) 状态下m 6 翼型的绕流进行了数值模 拟。本文采用的第一到第四层网格的网格数基本上是按8 ：1 递减的。第一到第四 层网格的单元数分别为：2 9 1 9 1 4 ，5 0 3 2 5 ，5 7 1 1 ，7 9 7 。图4 2 为三维m 6 机翼表 面各层网格的比较： 南京航空航天大学硕士学位论文 非结构多重网格法流场数值模拟 ( d ) 图4 2 采用的网格 f i g 4 2 t h ea d o p t e dm e s h e s 本章对单重算法及上文提到的v 型多重网格算法的收敛速度进行了比较结 果如图4 3 ，其中横坐标表示循环数，纵坐标表示残值误差的对数( 以1 0 为底取 对数) ，即收敛的量级。 oo o - 一10 0 l 九 i u 虬一二0 0 图4 3收敛率比较 f i g 4 3c o n v e r g e n c e r a t e 0 0 o 一 一 。一、 娟 j 南京航空航天大学硕士学位论文 5 0 0 0 05 0 一 一00 0 一 一05 0 0 0 i5 0 o 图4 4 等压线 f i g 4 4 p r e s s u r ec o i l o u r s 图4 5c 。分布( 4 4 呦 f i g 4 5 d i s t r i b u t i o no f t p ) 0 0 由图3 可见，跨音速( m = o 8 4 ，a ：3 0 6 。) 状态下m 6 翼型用多重算法只 需1 3 0 个循环便可以收敛到1 0 。5 量级，而单重算法需1 2 2 0 步以上。由于各层网 格数基本上是按照8 ：1 的比例递减，因此在粗网格上迭代花费的时间相对于第一 层密网格要少得多。因此用多重网格法比单重算法能大大的节约机时。另外，为 了验证多重网格法加速收敛的高效性，本文中除了多重网格法之外，没有采用其 它任何加速收敛的措施( 如残值光顺等) 。 非结构多重网格法流场数值模拟 图4 4 为跨音速( m 。= o 8 4 ，q = 3 0 6 。) 状态下m 6 翼型用多重网格算法得 到的机翼表面的等压线图。图4 5 给出了单重网格法和多重网格法计算的4 4 截 面处的c 。分布的比较，可见，结果完全重合。由于多重网格法中，在所有粗网格 层上进行的迭代只是为了更快地消除密网格上的低频误差，多重网格算法是否达 到收敛仍以密网格上的残值误差为标准，因此，事实上多重网格法计算的结果与 单重算法是一致的，只是收敛速度差别很大。 4 5 结论 本文设计了一种新的三维网格聚