淮洋，姚冰

（中國航空工業(yè)集團公司西安航空計算技術(shù)研究所，陜西西安，710065）

0 引言

在飛機概念設計中，涉及的幾何模型相對比較簡單，但需要大量的計算來評估飛機性能。引入CFD 數(shù)值方法，能夠在兼顧保真度和效率的同時，提供高精度的計算結(jié)果。網(wǎng)格生成是CFD 數(shù)值模擬的關鍵一環(huán)。

CFD 數(shù)值模擬飛機流場時，一般可采用結(jié)構(gòu)網(wǎng)格、非結(jié)構(gòu)網(wǎng)格[1-3]和混合網(wǎng)格[4-5]等策略。結(jié)構(gòu)網(wǎng)格策略面臨幾何拓撲限制和難以自動化生成的挑戰(zhàn)；而非結(jié)構(gòu)網(wǎng)格則面臨粘性網(wǎng)格生成的困難；混合網(wǎng)格策略彌補了純非結(jié)構(gòu)網(wǎng)格的不足?；旌暇W(wǎng)格將流場分為附面層區(qū)域和附面層外部區(qū)域。附面層內(nèi)部流動梯度大，各向異性三棱柱網(wǎng)格具有非結(jié)構(gòu)網(wǎng)格和結(jié)構(gòu)網(wǎng)格的雙重性質(zhì)，屬于半結(jié)構(gòu)網(wǎng)格類型，采用棱柱網(wǎng)格可以比較好地模擬飛機表面附面層效應；附面層外部流場變化相對緩慢，各向同性四面體網(wǎng)格無網(wǎng)格拓撲的限制，容易生成。因此，三棱柱/四面體混合網(wǎng)格策略能夠很好的兼顧網(wǎng)格自動化生成和飛機附流場的精確模擬。

現(xiàn)有的網(wǎng)格生成商業(yè)軟件（ICEM-CFD、PointWise、Gambit等）盡管功能強大，但都是通用型網(wǎng)格生成軟件，沒有專門針對概念設計中機身、機翼、發(fā)動機短艙、垂尾等特定部件的網(wǎng)格控制參數(shù)設定功能，無法保證網(wǎng)格的自動化快速生成。

本文發(fā)展了一種混合網(wǎng)格自動化生成技術(shù)：采用陣面推進法和節(jié)點光順方法自動化快速生成附面層內(nèi)部的棱柱網(wǎng)格；采用Denaulay、邊界恢復、網(wǎng)格細化等方法自動化快速生成附面層外部的四面體網(wǎng)格，用于CFD 數(shù)值模擬，為概念設計階段中飛機性能的快速評估提供支持。

1 三棱柱網(wǎng)格生成

三棱柱網(wǎng)格在附面層法線上有很高的網(wǎng)格分辨率和很好的網(wǎng)格正交性，又極大地減少了網(wǎng)格量，對于提高復雜外形的數(shù)值模擬精度和效率都是非常有利?；谌敲婢W(wǎng)格，采用陣面推進法和節(jié)點光順方法可實現(xiàn)附面層區(qū)域的三棱柱網(wǎng)格自動化生成，生成流程如圖1 所示。

圖1 三棱柱網(wǎng)格生成流程

1.1 節(jié)點推進矢量

節(jié)點推進矢量一般通過計算與該節(jié)點相連的所有面元法向矢量的平均值求得。

1.2 推進矢量光順

推進矢量光順公式為：

圖2 后緣角和凹角處推進矢量光順示意圖

1.3 推進步長

對物面上的點沿著推進矢量的方向進行推進，推進步長由下列公式進行控制：

式中，iδ為第i 層黏性網(wǎng)格的空間分布長度；1δ為預先給定的第1 層黏性網(wǎng)格的空間分布長度，它必須保證附面層內(nèi)部內(nèi)應包含若干節(jié)點；r代表預先給定的網(wǎng)格增長率因子。

2 四面體網(wǎng)格生成

四面體網(wǎng)格具有很大靈活性，可以對任意三維區(qū)域進行自動填充，便于實現(xiàn)網(wǎng)格自動生成。本文采用Delaunay 四面體剖分、邊界恢復、網(wǎng)格細化等算法實現(xiàn)自動化、快速地生成四面體網(wǎng)格，具體流程如圖 3 所示。

圖3 四面體網(wǎng)格生成流程

2.1 邊界恢復

網(wǎng)格模型中，邊界可用來區(qū)分模型內(nèi)外，或特殊物理性質(zhì)（如跡線、激波）的幾何位置。三維網(wǎng)格中，邊界通常定義為一組點、線和面。Delaunay 三角化只能保證給定點存在于四面體中，但不能保證給定線和面存在于四面體中。邊界恢復通過一系列局部操作修改初始四面體網(wǎng)格，使得最終網(wǎng)格中包含用戶指定的邊和面。邊界恢復包含兩步：先恢復邊界邊，再恢復邊界面。

在網(wǎng)格生成的邊界恢復以及優(yōu)化過程中會用到一些基本的變換運算，本節(jié)對這些基本運算做一簡單介紹。

運算（1）：在一條邊上插入一點，將以這一邊為邊的所有元素重新剖分，使元素增加一倍。

運算（2）：在一個面內(nèi)插入一點，將這個面相鄰的兩個元素重新剖分為6 個元素。

運算（3）：在一個元素內(nèi)部插入一點，將這個元素重新剖分為4 個元素。

運算（4）：在一個元素的兩個面上各插入一點，將這個元素重新剖分為5 個元素。這個運算相當于兩次進行運算(2)，但每次只對一個元素進行運算。

運算（5）：對一條邊周圍的所有元素進行運算，將這一條邊刪除，將留下的殼重新三角剖分。這個運算比較復雜，運算結(jié)果不唯一，甚至有可能得到非法剖分，因此采用時要十分注意。

運算（6）：對兩個相鄰三角形進行運算，刪除公共面，生成三個新元素。這個運算當兩個運算對象的組合為凸時可以進行，否則由于結(jié)果非法而不能進行。

運算（7）：對一點周圍的所有元素進行運算，將這一點移動一個位置，將所有元素重新組合生成相同數(shù)目的新元素。這一運算對點的移動位置有嚴格的要求，否則將會導致運算結(jié)果非法，運算失敗。

2.1.1 邊元邊界恢復

如果一條邊元AB 的丟失只與兩個元素有關。如圖 4 所示，其中AB 和兩個元素ACDE 和BCDE 相交。則只要將這兩個元素重新組合成三個元素ABCD，ABDE，ABEC 就可以恢復AB邊，這相當于執(zhí)行運算（6），由于三個元素組合成凸域，所以運算（6）可以執(zhí)行。一般情況下，AB 邊可能和多個四面體元素相交，可采用逐步計算的過程實現(xiàn)邊元邊界的恢復。

圖4 三維邊元破壞和恢復的簡單情況

2.1.2 面元邊界恢復

所有邊元恢復后，并不能保證所有面元均在網(wǎng)格中存在。圖 5 顯示了一個簡單的情況，其中三角形面元ABC 的三條邊均存在，但面元不存在，它被一條邊αβ 穿透。這種簡單情況下面元ABC 的恢復十分簡單，只要應用運算（5）于三個元素AαβB，BαβC，CαβA，得到ABCα 和ABCβ 二個元素就恢復了面元ABC。注意到αβ 和ABC 面有交點，三個元素組合為凸域，運算（5）可以執(zhí)行。

圖5 面元恢復的簡單情況

如果一個面元被多條邊穿透，則可以逐次應用運算（5）于每一條邊，使該條邊被刪除，最后恢復這個面元。

2.2 網(wǎng)格質(zhì)量

網(wǎng)格劃分中需要有質(zhì)量度量標準。一般網(wǎng)格線面夾角過大或者過小時，網(wǎng)格單元應該舍棄，因為這些網(wǎng)格單元降低數(shù)值計算的精度和效率。本文采用半徑邊緣比作為網(wǎng)格質(zhì)量度量標準。一個四面體t 有且只有一個外接球，令r 為四面體外接球的半徑，d 為四面體最短邊的長度，則這個四面體的半徑邊緣比ρ(t)為

式中: θmin 為四面體t 的最小面角。在判斷四面體網(wǎng)格質(zhì)量的應用中，半徑邊緣比要愈小愈好。大多數(shù)質(zhì)量差的四面體單元ρ(t)較大（如>2.0），除了sliver 型四面體。sliver 型四面體沒有很短的邊，但體積接近于0.

2.3 網(wǎng)格細化

在網(wǎng)格細化時，采用推進面插點法生成內(nèi)點。為了保證邊界附近生成較高質(zhì)量的網(wǎng)格，在選擇活躍的元素時，首先選擇具有邊界面的活躍元素，直到所有具有邊界面的活躍元素處理完，再處理其他活躍元素。網(wǎng)格密度的控制是由分布函數(shù)值控制網(wǎng)格外接球半徑，使小于給定的值。一般來說一個活躍元素的外接球半徑大于給定值，需要插入一點。插點過程中可見性的檢查是十分重要的。

如果整個計算過程中每一步產(chǎn)生的網(wǎng)格都是Delaunay三角剖分，則計算過程理論上不存在問題。由于恢復邊界邊元和面元的過程中可能會增加一些點，且造成部分網(wǎng)格不滿足外接球準則；計算過程的誤差也可能造成某些網(wǎng)格不嚴格滿足外接球準則。因此在生成新的內(nèi)點時進行一次濾波是必要的。即檢查新生成點附近是否已經(jīng)存在網(wǎng)格點。這一工作結(jié)合Delaunay 腔的生成過程進行是方便的，計算量也不很大。此外，在生成Delaunay 腔后還要檢查一下新點是否在腔內(nèi)，即是否屬于某一個元素內(nèi)。如果存在很接近新點的網(wǎng)格點，或Delaunay 腔不包含新點，則這個新點要舍去。

2.4 網(wǎng)格優(yōu)化

四面體網(wǎng)格在生成后，通常會包含一些網(wǎng)格質(zhì)量較差的單元，這對后續(xù)數(shù)值求解效、收斂性和精度有很大的影響。網(wǎng)格優(yōu)化能夠進一步提高網(wǎng)格質(zhì)量，因此非常重要。

網(wǎng)格優(yōu)化采用多種方式，如改變網(wǎng)格點位置，改變網(wǎng)格連接關系。最常用的是：網(wǎng)格點光順，邊/面交換，邊收縮，內(nèi)點插入。這些操作會根據(jù)具體情況相結(jié)合，迭代調(diào)整全局網(wǎng)格或者局部網(wǎng)格，使得網(wǎng)格單元的邊界縱橫比、最大最小二面角等參數(shù)不斷優(yōu)化，從而提高網(wǎng)格質(zhì)量。

3 算例測試

本文發(fā)展的混合網(wǎng)格自動化生成方法，在三角面網(wǎng)格的基礎上，通過簡單的參數(shù)設置，自動生成附面層內(nèi)部的棱柱網(wǎng)格和附面層外部的四面體網(wǎng)格。下面以噴氣式飛機模型為例說明混合網(wǎng)格生成的過程。

噴氣式飛機的物面網(wǎng)格為各向同性三角網(wǎng)格，在機身頭部、尾部，翼面前緣、后緣、梢部和根部等曲率較大區(qū)域進行了局部加密，物面網(wǎng)格單元規(guī)模約10 萬，如圖6 所示。

圖6 噴氣式飛機物面三角網(wǎng)格示意圖

三棱柱網(wǎng)格自動填充附面層內(nèi)部區(qū)域，通過設定第一層邊界高度和總附面層層數(shù)即可生成。圖7 給出了飛機展向?qū)ΨQ面上的三棱柱網(wǎng)格分布，無論在凹形還是凸型區(qū)域，均能生成較光滑的棱柱網(wǎng)格。

圖7 展向?qū)ΨQ面上三棱柱網(wǎng)格示意圖

該混合網(wǎng)格方法根據(jù)遠場中心和半徑參數(shù)，可自動化生成均勻分布的球形遠場邊界網(wǎng)格。四面體網(wǎng)格網(wǎng)格自動填充粘性邊界和遠場邊界之間的區(qū)域。從粘性附面層到遠場邊界的四面體網(wǎng)格可以按一定的網(wǎng)格增長率生成，特別的可以在模型外部建立一個加密區(qū)，使得加密區(qū)內(nèi)的網(wǎng)格尺寸小于設定值。圖8 給出了飛機展向?qū)ΨQ面上的四面體網(wǎng)格分布，從模型到遠場網(wǎng)格過渡較為光滑。

圖8 展向?qū)ΨQ面上四面體網(wǎng)格示意圖

表1 中給出了在相同的附面層初始高度、層數(shù)、增長比等控制參數(shù)下，本文方法和point wise 商業(yè)軟件的對比。基于該噴氣式飛機，由于幾何構(gòu)型比較簡單，本文方法的網(wǎng)格生成速度快于pointwise 軟件，在網(wǎng)格規(guī)模和網(wǎng)格質(zhì)量方面兩者相當。

表1 本文方法和商業(yè)軟件Pointwise 對比表

4 結(jié)論

本文發(fā)展的混合網(wǎng)格自動化生成技術(shù)能夠自動化生成附面層內(nèi)部的棱柱網(wǎng)格和附面層外部的四面體網(wǎng)格。采用噴氣式飛機三角面網(wǎng)格進行了驗證，在生成效率和網(wǎng)格質(zhì)量等方面與商業(yè)軟件水平相當。該混合網(wǎng)格自動化生成方法能夠在概念設計階段為飛機性能的快速評估提供技術(shù)支持。